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Abstract. The central objective of this paper is to develop reduced basis methods for parameter dependent 
transport dominated problems that are rigorously proven to exhibit rate-optimal performance when compared 
with the Kolmogorov n- widths of the solution sets. The central ingredient is the construction of computationally 
feasible "tight" surrogates which in turn are based on deriving a suitable well-conditioned variational formulation 
for the parameter dependent problem. The theoretical results are illustrated by numerical experiments for 
convection-diffusion and pure transport equations. In particular, the latter example sheds some light on the 
smoothness of the dependence of the solutions on the parameters. 



Over the past few years model order reduction has become an indispensable constituent of large scale design 
or optimization problems. In particular, the Reduced Basis Method (RBM) is perhaps by now one of the most 
important paradigms for highly complex frequent query problems involving parameter dependent PDEs, see e.g. 
[271131 [23]. Among other things, at least under certain circumstances, modeUng errors are rigorously controlled 
and can be upgraded if necessary. 

While the development of RBMs has been a very active area with impressive success stories in by now a 
variety of important application fields, it is fair to say that a theoretical underpinning of what one might call 
"near-optimal performance" - in a sense to be made precise later - is still confined to a relatively narrow problem 
class. The central purpose of this paper is therefore to extend the scope of problems for which RBMs can be 
developed and rigorously proven to perform in that near optimal sense. The focus of the present work is on 
performance in terms of the accuracy offered by the reduced model, roughly speaking, centering around the 
question how to ensure any certified target tolerance of the reduced model by a possibly small number of reduced 
basis functions, of course, always insisting on the standard offline- online division of the overall computational 
work. 

1.1. General Framework. Suppose that : X ^ Y' , fi E V, is a family of (linear) operators from a Hilbert 
space X onto the dual Y' of another Hilbert space Y, depending on parameters fi from a compact set V C MP. 
Under appropriate conditions on {B^J^g-p the solution set 



is a compact subset of X. In the context of frequent query problems, like steering a functional £{p{fi)) of the 
solution towards a target value, RBMs try to exploit the fact that A4 may be a very thin subset of X. In fact. 
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compactness of M means that the Kolmogorov n-widths 

(1.3) dn{M)x inf maxdist (TW, T/)x, 

dim V<n 

tend to zero as n — >■ c», where V is taken from the set of all 71-dimensional subspaces of X and 

maxdist (A^, X„)x sup inf 

The objective is then to construct (problem dependent) subspaces Xn C X of possibly small dimension n such 
that for a given target accuracy tol, say, 

(1.4) maxdist {M,Xn)x < tol 

is guaranteed to hold. In particular, this implies that for any p e and any bounded linear functional £ X' , 
a trivial estimate immediately gives \£{p) — i{Px,XriP)\ ^ ll^lltol (which could even be improved by duality 
arguments, see e.g. [4]), where Px,x„ is the X-orthogonal projection onto X„. 

Of course, a key question is how to practically construct spaces X„ warranting (1.4) for possibly small n. A 
common strategy of essentially all RBMs is the following. Given X„, find a surrogate R(ji,Xn), n G V, such 
that 

(1.5) Mti) - Px,xMl^)\\x < CRR{ti,X^) 

holds for some constant C/j independent of fi and n. Here it is crucial that the evaluation of R{p,,Xn) is 
sufficiently efficient so that the maximization of i?(/x, X„) over n G V is computationally feasible. Then perform 
the greedy algorithm GA based on this surrogate, as described in Algorithm [T] 



Algorithm 1 greedy algorithm 

1: function GA 

2: Set Xo {0}, n = 0, 

3: while argmax^gp RifJ--, Xn) > tol do 

4: 

fin+i :== argmaxi?(/^,X„), 

(1-6) Pn+l ■=P{f-ln+l), 

Xn+1 := span{X„,{p(/i„+i)}} = span{pi, . . . ,p„+i} 

5: end while 
6: end function 



We have ignored for the moment the fact that the snapshots p(/i„) can, of course, not be computed exactly 
but only approximately within some sufficiently large but finite dimensional "truth space" . 

To see whether such a greedy space search produces good reduced models one can compare them with the 



"best possible" spaces. Clearly, the n-width dn{M)x from (1.3) is a lower bound for the accuracy attainable 
by any RBM, i.e., 

(1.7) dn{M)x < crn{M)x sup \\p{fi) - Px,xMl^)\\x- 

Unfortunately, in general it seems to be impossible to compute the precise optimal subspaces for which the 
71-width is attained. Nevertheless, the closer an{A4)x is to dn{M)x the better the choice of Xn- 

To see what can be achieved in this regard, recall from [H [3] that even when _R*(/i,A"„) := ||p(/i) — 
Px,x„p{^J')\\x is the ideal surrogate, in a direct comparison cr„(A^)x < Kndn{M)x the constant Kn can 
be as large as 2". Nevertheless, the following more favorable results in terms of convergence rates hold for 



surrogates that are tight, i.e., if in addition to the upper bound ( 1.5 ) it uniformly sandwiches the exact distance. 
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Definition 1.1. We call the surrogate R(fj,,Xn), tight if there exist positive constants cr^Cr such that 
(1-8) Cfli?(Ai,X„) < \\p{^i) ~ Px.xMWx < CRRifi,X„), 

uniformly in fi V . Moreover, we call 



(1.9) k{R) := inf {Cr/cr : cr, Cr satisfiy ([Ts]) for all^i£V,n£ N}, 

the condition of the surrogate R. 

Remark 1.2. As already observed in whenever the surrogate is tight, i.e. (|1.8[) ) holds, then the snapshots 



Pn — pifJ-n) from (1.6 1 satisfy the weak greedy condition 



(1.10) \\pn~ Px,x„Pn\\x >K{R)-^msixdist{M,Xn)x, neN, 



where k{R) is given by (1.9 I 



The following statements are then readily derived from the results in pQ [TT] . 



Theorem 1.3. Assume that the spaces Xn are obtained through a greedy algorithm GA, (1.6) based on tight 



surrogates. Then, if dn{M)x = 0{n ") , for some a > Q or if dn{A4) x ~ 0{e ™ ), for some c, a > 0, one has 

(1.11) maxdist (7W,X„)x = 0(n""), maxdist (7W,X„)x 0(e^™°), n -J> oo, 

respectively, where the constants depend on a,c, and k{R) with exact specification given in [UlTTj. Moreover, 
these bounds remain valid up to the tolerance tol* when all computations are carried out within this accuracy. 

We call an RBM rate-optimal if the generated spaces Xn satisfy "Kolmogorov optimal" bounds of the type 

There are two important points to be drawn from these results that guide the subsequent developments. The 
first one is: although dispensing with the (infeasible) ideal surrogate i?*(/x,X„) := \\p{fJ^) — Px,x„p{tJ-)\\x, a 
tight surrogate still ensures that the accuracy provided by the reduced bases is in terms of rates still essentially 
as good as that of the "Kolmogorov-best" subspaces. The second point is quantitative. It is absolutely vital to 
make sure that the condition k{R) stays as small as possible. In fact, a look at the dependence of the constants 



in (1.11) on k{R) (see [U [TT]) reveals that the closer k{R) is kept to one, the better is the accuracy of the 
reduced spaces, in comparison with the best spaces, already for a small reduced dimension, which is at the heart 
of model reduction. 

Hence, the central objective of this paper is to develop a rigorous conceptual framework to obtain practically 
feasible tight surrogates whose condition k{R) < Cr/cr is as close to one as possible, in particular, for problem 
classes for which this is currently not known. 

1.2. Objectives and Layout. To provide an orientation for subsequent developments the corresponding ideal 
scenario and the corresponding basic mechanisms are briefly recalled in Section [2. 2 [ It is by and large confined 
to problems that are uniformly elliptic with respect to the parameters. The perhaps next best understood 
case is the reduction of a parabolic problem to a sequence of elliptic problems [T71 [Ul [TB], where however. 



the lower bound in (1.8) - and hence tightness - seems to be missing. This has been recently significantly 
improved in |24j using a space-time variational formulation. Moreover, important progress has been made in 
[13l[T4l[26] developing RBMs for specific saddle point - hence indefinite - problems such as the Stokes system. In 
particular, in the present terminology stability and, as a consequence, tight surrogates are obtained by enriching 
the velocity spaces by supremizers. More precisely, there are two approaches. For standard affine parameter 
dependence of the involved bilinear forms one can determine a priori an enrichment, depending on the number 
of terms in the bilinear forms, that ensures that the infinite dimensional inf-sup-constant is preserved, see in 
[13l|26]. Since the number of these supremizers is possibly quite large, as an alternative, it is proposed in [14] 
to adaptively add supremizers until a desired inf-sup-stability is reached. It is observed experimentally that in 
the tested examples this adaptive enrichment results in an overall much smaller number of stabilizing functions 
although the actual guaranteed termination of such a procedure has apparently not been discussed. Although 
termination in the context treated in [131 114L I26| is apparent, we shall encounter situations where this is no 
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longer the case. Nevertheless, relating also the stabilizing enrichments to greedy approximations allows us to 



treat this case as well, see Section 4.4 



Although the present paper addresses a rather different problem class the treatment of saddle point problems 
turns out to be an important point of contact. In fact, the stabilizing enrichment of the reduced velocity spaces 
by adding supremizers can be viewed as a special instance of the interior loop of what we call here double 
greedy schemes^ presented first at a workshop in Paris, 2011 [H]. The central objective of this paper is in fact 
to develop rate-optimal RBMs - viz. identify well-conditioned tight surrogates - for a much wider scope of 
problems, in particular, to those that are at present notoriously not covered by current RBM theory, namely 
transport dominated problems. Two model problems are formulated in Sections [2.3.H|2.3^ exhibiting increasing 
levels of obstructions. We emphasize though that the general methodology presented below is not restricted 
to those problems at all. A key role in this context is played by deriving stable variational formulations for 
such problems that are necessarily of Petrov-Galerkin type. The main features of this approach, being valid for 
a wide range of problems including indefinite, unsymmetric and singularly perturbed problems, are shown in 
Section [3] They can be summarized as follows: 

(i) Tight a-posteriori bounds for the truth spaces as well as reduced spaces warrant certification. In 
particular, truth and reduced spaces can be upgraded without discarding prior computations, see the 
robustness results in [T]. 

(ii) While remaining feasible in the sense of an onlinc/offiine decomposition through a built in stabilization 
loop, the scheme automatically gives rise to stability constants that can, in principle, be made arbitrarily 
close to one, see also (5.7) and Section [6j 

(iii) Viewing time as an additional "spatial" variable, the results can be applied to time dependent problems 
through corresponding space-time discretizations, which is one reason to focus on transport problems, 
see [7]. 

In summary, the particular variational formulations presented in Section [3] combined with certain stabilization 
techniques optimally inherits the analytic structure of the underlying infinite dimensional problem to the reduced 
model. 

Section [5] is then devoted to the algorithmic development and analysis of a double greedy scheme giving rise 
to rate-optimal RBMs. 

The theoretical findings are then applied in Section [6] to the two model problems concerning convection- 
diffusion and pure transport equations. First numerical experiments quantify the results and highlight several 
particular obstructions. 

In Section [t] we apply the (slightly modified) scheme to other types of saddle point problems not necessarily 
stemming from the generation of well-conditioned variational formulations. As a simple consequence we obtain 
rate-optimality also for the problems considered in [131 Ull ES] . 

To simplify the exposition we write a < 6 to express that a is bounded by some constant multiple of h 
independent on any parameters a, b may depend on. Likewise a ~ & means a < b and b < a. 



2. Conceptual Preview 

2.1. Feasibility. In all subsequent developments we will be dealing exclusively with affine parameter depen- 
dence, see e.g. .23- Under this assumption we insist on the usual division of the computational work into an 
offline and online mode. Solving a problem in the full space X, which is typically computationally very intense, 
happens only in offiine mode where it is understood that actual computations take place in some sufficiently 
large but finite dimensional subspace of X which is commonly referred to as the "truth space" . Typically 
X^ is chosen so as to guarantee 

(2.1) sup inf \\p-v\\x < tol*, 

for some tolerance tol* that is sufficiently small for the application at hand. The subscript M is sometimes 
supressed when there is no risk of confusion. The greedy search for the reduced basis functions falls therefore 
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into the offline mode. Tliis requires evaluating the surrogate for a sufficiently large training set of parameters 
which for simplicity we also denote by V. In what follows, we call the surrogate feasible if each evaluation of 
the surrogate requires solving only a problem in the small current reduced space X„. We sometimes say then 
that the offline mode is (computationally offline) feasible. 

Likewise, the online evaluation is called feasible if each reduced basis approximation of some p{p) requires 
solving only a "small" problem of dimension n in the reduced space Xn ■ In this mode solving a "large" problem 
in is prohibited. 

Note that a feasible surrogate is not allowed to explicitly contain the true solution (in the truth space). 
This is why one is essentially forced to resort to residuals to estimate the true error, which in turn requires a 
tight error-residual relation. 

2.2. The Ideal Setting. The following "ideal setting" reveals the basic mechanisms leading to residual based 
tight surrogates. 

To this end, let 6^(-, ■) : X x X ^ R he a. symmetric uniformly X-elliptic bilinear form and £ E X' , i.e. 
(2.2) Ca\\q\\jc<b,,{q,q), bf,{p, q) < Ca\\p\\x\\q\\x , P, q e X , fi € V , 

holds uniformly in ^ E V. For compact V one obtains a compact solution set C X for: given i E X' , find 
j3(/i) E X such that 

b^,{p{n),q) ^ {i,q), qEX. 
There are two key properties that ensure rate-optimality in this setting: 



(MP) Mapping property: the operator _B^, defined by {Bfj_p,q) = 6^(p, g), p,q E X, is for each ^ e 7^ an 
isomorphism from X onto X' , i.e. 

(2.3) \\p\\x ^ WBfj^pWx', uniformly m eV. 

In other words, errors measured in the "energy norm" |j • |jx are equivalent to residuals in the dual norm || • \\x'- 
(BAP) Best Approximation Property: The Galerkin projection to the current reduced space, which can 
be done in online mode, produces, up to constants, a best approximation with respect to the X-norm. 
In fact, denoting by x„ the Galer kin-projector onto X„ defined by 

bt,{p{f^),q) = b^,{U^,xMtJ-),'l), 1 e Xn, 
combining Cea's Lemma with MP provides for p„(/i) :— Wf^^x^Pil^) 

fo A\ W I \ V ^ Ml {^,q) -bf,{pn{fi),q) , . 

(2.4) \\p{ii) - Px.xMnWx ^ sup :=i?(/i,X„). 

qex \\q\\x 



Thus, MP and BAP imply that the residual based surrogate, defined by (2.4), is tight, while the computation 
of pip.) is completely avoided but traded against the cheap computation of the Galerkin projection in X„. 
However, the condition k{R) of the surrogate (see (|1.9[)) depends on the condition number Kx.x'iBfj} < Ca/ca 



(see (2.2)) of the operator i?^, which should therefore be of moderate size. 



Finally, feasibility of the surrogate in (2.4) is well known (see e.g. [57]) to be ensured when the parameter 
dependence of 6^(-, •) is affine, see (3.42) below. 



2.3. Two Model Problems. As soon as one leaves the elliptic setting MP, BAP, tightness of residual 
based surrogates, are no longer for free. In particular, so far well-conditioned tight surrogates do not seem 
to be available yet for many unsymmetric PDEs like convection dominated or pure transport problems. We 
shall discuss two model problems that bring out several principal obstructions. The first example concerns 
convection-diffusion equations for which, in principle, classical variational formulations are available. The 
second example concerns pure transport for which a "natural" variational formulation is less obvious and for 
which the parameter dependence of the solutions turns out to be less regular. Perhaps more importantly, the 
two examples represent two different scenarios regarding the spaces associated with the bilinear form b^{-,-), an 
issue that has apparently not been addressed in the REM context. 
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2.3.1. Convection-Diffusion Equations. As a first example we consider the linear convection- diffusion equation 

(2.5) -div(eVp(a;)) + 6(^)- V]5(a;)+cp(x) = /(a;), in^, p = on dn, 

where for simplicity we assume for now that only the convection depends on a parameter while e could 
be arbitrarily small. We could as well include the viscosity and the reaction term varying in suitable regimes. 
Its classical weak formulation is 

(2.6) 6^(p,<z):=e(Vp,Vg) + (%)-Vp,(?) + (cp,(z)-(/,g), q^X^Hlin). 
It is weU known that h^if) G W^^'^iVlY, c e L^{n), e V, such that 

(2.7) - ^divfe(^) +c > 0, 



implies well-posedness of (2.61 in the sense that the induced operator : Hq (£7) — > {Hq (fJ))' is an isomorphism, 



i.e., there exists for each ^ e P a unique solution p(/i) to (2.61 in ffQ(r2). However, although (2.2) is still valid, 
the condition number K//i(n),H|J(n)(^A') behaves like the Peclet number |5(/Li)|oo/e and hence is inacceptably large 
for strongly dominating convection. As a consequence, in this case the condition k{R) of the corresponding 



surrogate (2.4 1 based on the H ^(ri)-residual grows with the Peclet number. Hence, although such a surrogate 



is theoretically tight, as long as e > eq where eo > is fixed, the condition k{R) (see (1.9)) is so large, that, due 



to the constants in (1.11), one can expect essentially no control of the quality of the reduced spaces for very 
small eo and moderate n. 

Therefore, we are mainly interested here in a robust treatment of arbitrarily large Peclet numbers |6(/i)|oo/e 
which to our knowledge is currently not well covered by RBM methodology. 

Unfortunately, an easy cure based on the standard (mesh-dependent) stabilization methods such as SUPG 
(see e.g. [25] for a survey) does not give rise to an error-residual relation that stays independent of the Peclet 
number |6|oo/e either. 

Instead we pursue here a different line based on stabilizing the problem on the infinite dimensional level 
which, in particular, involves unsymmetric variational formulations, i.e., 6^(-, •) is viewed as a bilinear form on 
a pair of (possibly) different and parameter dependent Hilbert spaces X^,y^, ii CzV. 

2.3.2. Linear Transport Equations. In some sense the situation is even aggrivated when the diffusion vanishes 
completely as in pure parametric transport equations forming the core ingredient of Boltzmann equations and 
related kinetic models as well as kinetic formulations of conservation laws. Already the simplest version of a 
(stationary) linear transport equation 

(2.8) f,-Vp + cp^f, inr!=[0,l]^ p ^ Pb, onT^i^i), 

will be seen to represent the "worst scenario" from the RBM perspective, where, denoting by n{x) the outward 
normal at the point x, 

r_(Ai) := {xedfl: n{x) ■ fi<0}, 

is the inflow boundary for the given convection vector /x. An example of a parameter domain would be the 
sphere S'^~^ appearing in radiative transfer models, see [HI [28]. It will be seen that the two examples differ 
in a subtle but essential way, in particular, regarding smoothness of the dependence of the solutions on the 
parameter. 



A possible variational formulation of (2.8) can be found in [12j. In order to eventually apply the n- width 



benchmark, it is preferable to measure all parameter dependent solutions in a single reference norm. Therefore 



we employ here a slightly different variational formulation from [7 : multiplying (2.8) by a test function and 
integrating by parts, yields 

(p, -^i-Vq + cq)+ j n- ^pq ^ (f^q) - n- fipq. 
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If we now take test functions q that vanish on dfl \ r_ the boundary integral on the left hand side is zero. 
Furthermore, we may replace the function p in the boundary integral on the right hand side by the boundary 
condition pb so that we obtain 



(2.9) b^{p,q) := {p,-^i-Vq + cq) = {f,q) - J n ■ fipbq. 

For this variational formulation it is natural to define the function spaces 

(2.10) clos||.||,^^ {q e C°°{n) : q\en\r^ = O} , := L2(f^) 
endowed with the norms 

(2.11) My, WBZqh,, \\p\\x, := ||p||l,. 



It is shown in that the operator induced by 6^ (.,.), is an isomorphism : — > Y"^ so that (2.11) 
indeed defines a norm. 

Remark 2.1. Notice that the spaces differ even as sets for different ^. Moreover, in contrast to the previous 
example we must have Xf^ ^ here. 

3. Robust Error- Residual Mappings for Unsymmetric Problems 

3.1. The Basic Principle for MP. In the following, we consider general bilinear forms .) : x Y"^ ^• 
M for possibly parameter dependent Hilbert spaces and giving rise to what one may call an infinite 
dimensional Petrov-Galerkin formulation where the trial space X^ generally differs from the test space 1"^. 
Thus, the operator given by {B^q, v) = h^{q, v), q £ X^, v £Y^,\s now viewed as a mapping from X^j_ to 
Y^. In accordance with the preceding examples, we shall assume that this operator is actually an isomorphism, 
i.e., the operator equation 

(3.1) B^p{^i) = /, 

has for any / £ Y^' a unique solution in X^. However, may possibly have an unacceptably large p,- 
dependent condition which can be quantified with the aid of Babuska's Theorem: if there exist constants 
< P{p),Cb{fi) < oo such that 

o^ ■ f bf,{q,v) \bf,{q,v)\ ^ 

(3.2) mf sup 11 11 |. 11 > /3(m), sup sup , < Cb{n) 
q&x,^^Y, M\xJv\\y, qex,veY, llgllxjiwil y,, 

and for every v G there exists a. q £ such that b^{q,v) ^ 0, then one has nx,.Y,{B^) < Cb{fi)/ P{fJ,). 

As in the case of convection dominated convection-diffusion equations kx ,v', (S^) could be very large, severely 
degrading a greedy selection of snapshots in a REM. The goal of this subsection is to describe how to "stabilize" 



the problem on the infinite dimensional level which could be viewed as preconditioning (3.1). The underlying 
basic principle has been used before in several works for different purposes, see e.g. [51 HOI 1121 122]. Here our 
main orientation is taken from [51 [71 [32] . We briefly rehash the essential facts in order to bring in an additional 
new element, namely the interrelation of Petrov-Galerin schemes and associated saddle point problems, which 
plays an essential role for eventually constructing well-conditioned tight surrogates. 

We begin with collecting a few useful preliminaries. It will be usefull to identify for a given q £ the 
supremizer Vq for which sup^^y^^ bf^{q,v)/\\v\\Y, is attained, see e.g. (23i il3j . 

Remark 3.1. For every q G X^ the optimal test function is given by 

o\ bf,{q,v) 1 

(3.3) Vq := argmax-j— n = Ry B^q, 

veY, \\v\\y, " 

where Ry, :Y^ ^ Y^ is the Riesz-map defined by 

(3.4) {Ry^v,w) = iv,w)y^, v,w&Y^, \\-\\1^{;-)y,, 
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Hence, in particular, one has 



mi T—r = ini sup 



hWx^ <l'^X^^veY,A\q\\x^\v\\Y,, 



For convenience we recall the simple argument. Written in variational form, the supremizer is defined by 
(vq,w)Y^ = bfj.{q,w) for all w G Y^, which yields 



bt,{q,v) {vq,v)Y^ 
veY^ My^ veY^ \\v\\y^ 



which readily confirms the claim. 



Although for most of the following considerations the dependence of the involved bilinear forms on the pa- 
rameter /i € 7^ is irrelevant it will be convenient for later purposes to retain the parameter dependence in the 
notation. 



Renormation: The possible ill-conditioning refiected by a very large kx,y{B) < Ct/P in (3.2 ) can be remedied 
by properly modifying one of the two norms || • ||y, or || • \\x while keeping the other one fixed. Specifically, we 
wish to choose an equivalent but possibly different norm || • ||^ for X so that ideally Cf,(/i) = /3(/i) = 1, which 
then means that 

(3.5) - llZ-Sp^lly;, qeX,,^ieV. 

In this event the residual of a best approximation would be an ideal surrogate even sparing one the computation 
of stability constants for the error certification. 

Our starting point is exactly this latter ideal error-residual relation. Specifically, given || • \\y^, we endow now 
Xf^i^ with a new norm |1 • |1 j,^ j defined by 

(3.6) Ibb, ^= sup - - \\Ry'B^p\\Y,, p e x„ ^^ev. 

veY^ \\v\\y^ ^ " 

Note that this is indeed a well-defined norm because Ry^B^ : 1^ is an isomorphism, hence injective, 

and that because of WRy^^B^qWy^ = {B'^Ry^B^q, q) , the corresponding Riesz map Rj^^ : X^, — )• X'^ is given by 

Rj^ :— B*^Ry^ B^. In addition, this shows that the X^-norm is equivalent to the original norm, i.e., there are 
Cm , C'j\/ > such that 

(3.7) CM\\q\\x,<\\q\\x<CM\\q\\x,. qeX^,iieV. 



Note also that (2.11) is a special case of (3.6), where 

(3.8) Ry, - b^b;. 

and thus II . 11^^ = \\{B^B;)-^B^ • l|y, = II • Ik.. 

Remark 3.2. For the \\-\\x norm one has optimal continuity and stability constants Cb{fJ.) = /3(m) = 1, /i G 'P, 
i.e. 

(3.9) sup sup I — = inf sup — = 1. 

qex^vGY^ WnyJqWx^ lex^veH^ ll^llyjkllx,. 

Hence, k-^ yi (^m) ^ '■^ isometry for these norms, which is the desired robust - in fact optimal - 

error-residual relation (3.5) MP. 
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Proof. The first relation follows from 

Mq,v)\ = \{Ry'^B^q,RY^v)\<\\Ry'^B^q\\Yj\RY,v\\Y^ 



On the other hand, note that for any q G its supremizer Vq := Ry^B^q G gives by (3.6|, (3.4), 

btj,{q,Vci) = {B^,q,RyBf,q) = 



(3.10) 



which completes the proof. 



li and lli^^lly^ = 



inf sup 

qex^y^Y \\q\ 



l\\ so that 

^ . , {B^^qMylB.q) 
> mi ,, = 1, 



qex^ 



□ 



3.2. Petrov-Galerkin and Saddle Point Problems. The validity of BAP is no longer automatic for un- 
symmetric formulations. In principle, it can be approached through contriving suitable Petrov Galerkin dis- 
cretizations. A central issue in this section is to relate such Petrov-Galerkin schemes to equivalent saddle-point 
problems. In particular, this avoids the explicit computation of the respective test spaces which could be 
parameter dependent. 

To this end let W C A^ be a "generic" trial space which will play several different roles. It may stand for 
the full infinite dimensional space, or for the truth space, or eventually for the reduced space. Notice first that 
the best approximation Pwip) € W for p{fi) = B~^ f is 

(3.11) Piv(m) := argmin||p(^) - (?||^ = argmin ||/ - y/ , 

qGW ^ qeW 

which is therefore given by the normal equation: find pwip) G W such that 

(3.12) {f-B^pw{tL),B^q)y,=0, qeW. 

What keeps us from using this as the basis for a variational discretization, is the fact that the y^-scalar product 
is usually hard to evaluate numerically. Noting that Ry^ — Ry^ the last equation is equivalent to 

Ry' if - B^^pwiii)) , B^q) =0, qeW. 



(3.13) 

Introducing the auxiliary variable u(/i) :— Ry' (/ — B^pw{^)), or rather 
(3.14) {RyAi'),v) = {f - B^pw{ti),v) , v€Y^, 



in weak form, the relation (|3.13|) and hence (|3.11|) can be equivalently written as 
(3.15) 



{RYu{^i),v)+bf,{pw{fJ.),v) = {f,v), veYf,, 



b^{q,u{fj,)) 



= 0, 



qew, 



which now just involves standard L2-inner products. Of course, in particular for W — 



(3.16) 



{Ry^u{h),v) +b^{p{^i),v) 
bf,iq,u{n)) 



0, 



is equivalent to the original problem (3.1), which now takes the form of a saddle point problem. Bijectivity of 
Bf^ readily shows that 

(3.17) u{p) = 0, neV. 



Hence, the solution manifold of the saddle point problem (3.16) in [J^^-p A^ x can be identified according 
to 



(3.18) 



M=Mxx {0}, Mx ■■= {p(/i) : p{fi) solves (3.1)}, 



as it should, with the one for the original problem (1.1). 
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Remark 3.3. Even when the spaces Yf^ differ as sets when jj, varies, as e.g. in (2.11 1 for the transport equation, 
the solution manifold is still compact as long as the norms \\ ■ |jjsf^ are all equivalent to a reference norm. Hence, 
the greedy errors are guaranteed to tend to zero and the n-widths benchmark is applicable. The issue of parameter 
dependence of the involved spaces will be taken up in Section \3.3\ again. 



Now given a finite dimensional subspace W, we cannot treat (3.15) yet, since we cannot test by all u e Y^. 
The following interpretation of this idealized situation is immediate from the normal equation (3.13). 



Remark 3.4. The problem (3.151 is equivalent to the Petrov-Galerkin scheme: findpwip) such that 

(3.19) bfj,{pw{l^),v) ^ {f,v), veYw, 
where 

(3.20) Yw 



is the optimal test space associated with W and pwi^J■) is the best X ^-approximation to p{^i) in X^^ 



Since (3.19) is practically infeasible a natural strategy is to replace by a sufficiently large finite dimensional 
subspace V GY^ that inherits "sufficient" stability. The following observation, which plays a crucial role in what 
follows, explains the interrelation between a practically feasible version of (3.19) and a fully finite dimensional 
version of (3.15). 

Proposition 3.5. The solution component pw^v (fJ-) of the saddle point problem 

{RY^uv,w{fJ'),v) +bf,{pw,v{fJ.),v) = veV, 
b^_i{q,uv,w{lJ')) = 0, q^W. 



(3.21) 



solves the Petrov-Galerkin problem (3.19) with the optimal test space Yw replaced by Yw = Py^.v{Ry^ Bf,{W)) 
where Py^^.v denotes the Yf^-orthogonal projection. 



Proof. For any q S W, consider Vq := Py^yiRy Bfj.l) G ^ and note that, by the first equation (3.21 1, 



bf,ipwyilJ-),Vq) 



{B^PW,v{lj),Vq) = {f,Vq) - {RY^Uv,W,Vq)- 



Since 



{Ry^Uv,W,Vq) = {uv,w{lJ'),Vq)y^^ = {uv,w{lJ.),Ry^Bf,q)y^ ^ bf,{q, Uy^w (l^)) = 0, 

where we have used the second equation in ( |3.21 ). 



□ 



Clearly, the larger V the closer Yw is to Yw so that the choice of V can be viewed as a stabilization. To 
quantify this observation, we call V is S-proximal for W if 

(3.22) 11(7 - PY^^y)Ry'^B,q\\Y^ < S\\Ry'^B^q\\Y^, q € W, 

holds for some fixed < S < I. see [7ll32]. 



Proposition 3.6. Assume that for given W x V C X^ x Y^ the test space V is 5 -proximal for W, i.e. (3.22) 
is satisfied. Then, one has 

1 



(3.23) 

and 
(3.24) 

Moreover, one has 
(3.25) 



\\p{l^)~Pw,v{l^)\\ 



< 



I - 6 qew 



lb(M) -mv(M)ll 



I"(m) - wy(/i)||y < inf ||p(/i) - q\\ 

1 — qew 

bp,{q,v) 



inf sup 
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Proof. Let pwip-) denote the best v\"^-approxiination to the exact sohition of (3.161. Then, for any q Cz W 
one has, on account of Remark |3.2[ 



{PwifJ') -Pw,v{f^),q)x^ = (p(/^) -Pwy{jj),q)x^ = {B^,{p{fJ,) -pw,v{p)),B^q)Y^ 

= {B,,{p{^^) - Pw,v{f^)),RY^^{B^,q)) ^ h^{p(p)-pw,v{ij),RYl{Bt,q)) 

= bf,{p{fi) - pw,v{fJ-), (I - PY^,,v)RY^{B^q)), 

where we have used Petrov-Galerkin orthogonahty, asserted by Proposition |3.5| in the last step. By duality. 
Remark 3.2 (3.5|, respectively (3.6), and ( 3.22[ ), we conclude that 

„ , , , ,„ b,ip{p)~PW,v{p)Al-PY,y)RY'AB,q)) 

\\Pw[N~Pw,v[mx, 



= sup 
qew 



< 



PifJ-) ~ Pw,v{^J')\\xJ\\^Y^iB^' 



S\\pM -Pw,vM\\ 



from which (3.231 follows by triangle inequality. 

Next recall from (3.171 that, in view of the first relation in ( |3.21[ ), 

\W{^J) - Uv,w{^J)\\Y,^ = \\uvm{P')\\y^ ^ {uv^w{fJ'),uv,w{fJ'))Y^ = {.f - B^pw,v{^i),uvM{p)) 

< \\f - B^pwy{tJ)\\Y'\\uv,w(.l^)\\Y^ = \\p{^^) - Pwy{^J.)\\x\\u{^x) - uv^w{^^)\\Y^, 



which together with (3.231 confirms (3.24) 



Finally, the inf-sup estimate ( 3.25 1 is an immediate consequence of the more general Proposition 3.7 below. □ 



We shall use the saddle point formulations to contrive rate-optimal RBMs, namely, on one hand, for comput- 
ing truth snapshots in = Xj\f with a suitable 5-proximal test space Ya/i a-nd on the other hand, for computing 
Galerkin projections in reduced spaces W = Xn again with an associated (5-proximal test space y„, whose 
construction will be discussed in the next section. 

To put this into proper perspective, given any W C X^, the condition ( |3.22 1 on a y C implies the 
best approximation property BAP for the Galerkin solution component pw.vifJ-) of (3.21) with a constant 
that becomes the closer to one the smaller the relative error becomes in ( 3.22[ ). Moreover, (3.24) says that the 
accuracy of the second "auxiliary" component uv,wil^) is automatically completely governed by the accuracy 



of the first component pw.vifJ^)- Finally, (3.22) implies inf-sup stability of (3.21). It will be shown below (for 
later purposes in a little more generality) that conversely inf-sup stability (3.25) implies (5-proximality. In fact, 
since the bilinear form a^(u,i«) := {Ry^v,w) = {v,w)y^ is trivially y^-elliptic with coercivity and continuity 
constants Ca{fJ,) = Ca{fJ.) = 1 (see (2.2)), we could have derived the best approximation property BAP (3.24) 



directly from a uniform inf-sup condition from standard facts about general saddle point problems, see e.g. 
We have presented the relatively short self-contained derivation in order to identify the precise constants and to 



bring out the particular role of the (5-proximality condition (3.22). As we shall show later both conditions (3.22) 



and (3.25) can be used algorithmically to ensure stability of the saddle point problem and hence tightness of 



corresponding residual based surrogates. 



The above discussion draws essentially on the use of the particular norm (3.6). In the context of classical 



saddle point problems, such as the Stokes system, it is more convenient to work with the "original" X^-norm 
related to II • l|x by (3.7). The following proposition clarifies the announced interrelation between an inf-sup 



condition and (5-proximality. 

Proposition 3.7. As before assume that y C ^'^'^ W C X^, < S < I, ana 
conditions: 

(1) 
(3.26) 



A > 0. Consider the two 



11(7 - PY,xJRY'B^q\\Y,. < S\\R-^B,q\\Y^, ^qeW, 
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(2) 
(3.27) 



inf sup 



qew^^v \\v\\Yjq\\x,, 



> A. 



Then (i) implies (ii) with constant A = cm Vl — 5^ ■ Conversely, (ii) implies (i) with constant 5 — ^1 — C^X^ 
i.e., A = C^/ \/l — 5'^ , where cm o.nd Cm o-'^e the constants form the norm equivalence (3.7). 



Note that when cm = Cm = 1, e.g. in case we use the X^-norm for X^, both stabihty conditions are 
equivalent. 

Proof. We reformulate (i) and (ii) in terms of equivalent conditions that can be more easily compared. First, 
squaring (3.261 and using that -Py^^.y is the y^-orthogonal projector, we obtain 



\PY,yRylB^q\\\<5^\\RylB^q\\l^, 
which is equivalent to 

^l^\\RylB^q\\Y, < \\PY,.yRy'^B^q\\Y,, 

By the definition (3.6) of the graph norm |1 ■ ||x j this is equivalent to 
(3.28) WPy^yRy'^B^qWy^ > ^l^Wq^ 



Vg e W, 

\/q e W. 
Vq e W. 



Next, to reformulate (3.27). Obviously the inf-sup condition is equivalent to 



(3.29) 



sup 

vev \\v\\Y^ 



>M\q\ 



x„ 



Vq e W. 



From (3.3) in Remark 



which yields 



3.1 



(3.30) sup 



bti{q,v) 



we know that the left hand side is maximized by the function v — PY^.vRy^B^q 

[Py^yRylB^q, B^q) {Ry^Py^yRy^B^q, Ry^B^q) 



vev miY, 



IPy^yRy'^B.qWy^ 



\PY^yRY^Bf_,q\\y^ 



= \\PY,^vRY'B^q\\Y^ 



Substituting the right hand side in the left hand side of the condition (|3.29|), yields 

Vg e W. 

|x„ and 



\PY^yRy'B.q\\Y>X\\q\\x,,, 



We see that this is condition is identical to (3.28) up to an equivalence of the 
proves the assertion. 



norms, which 
□ 



In summary, given trial space W C X^, a suitable V CY^ such that the Galerkin problem (3.21 ) has the best 



approximation property BAP, thereby warranting tight residual based surrogates, can be obtained by realizing 

b,,{q,v) 



(3.31) 



qew„^v IHYjqWx, 



where (3 := min^gp /3(/x) > 0, see (3.2) 



3.3. Parameter Dependence, Truth Spaces, and Feasibility. Before applying the above findings to the 
construction of well-conditioned tight surrogates, we need to be a bit more precise about the parameter depen- 
dence in order to distinguish eventually several relevant scenarios. Notice that the spaces Y^,X^ are allowed 
to depend on fi £ V in a way that they even differ as sets and no parameter independent reference norm may 
exist, see Remark |2.1| Let 

(3.32) r := fl X :^ f]X^, 

where the intersection is understood in the sense of sets. It is clear that Y and X are linear spaces. Although 
in general, we do not insist though, that Y and X are endowed with norms that are equivalent to all || . ||y and 
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II • \\x^, respectively. However, we do assume in what follows that Y,X are dense in Y^, X^, respectively, for all 
^ € V. Moreover, on account of the compactness of T-", wc can always define (possibly stronger) norms 

(3.33) 



\\v\\y := sup ||i;||y,,, 



mix ■■= sup \\q\\x^ 



for Y,X, respectively. Moreover, since Y = H^gp is assumed to be dense in Y^, for the inf-sup condition 
(3.2) it suffices to take for V — Y^ the supremum over Y instead of F^, i.e., there exist subspaces F C F for 
which the discrete inf-sup condition (3.27) holds uniformly in the parameter /i. 

Of course, this setting covers, in particular, the special situation - usually considered in the RBM context 
- that all the spaces Y^,X^, fi €V, agree as sets, respectively, and where the respective norms are uniformly 
equivalent, i.e., there exist constants < Co, Co < oo such that 



(3.34) 

and 

(3.35) 



collt'lly < My, <Co\\v\\y, ^leV,veY, 
CohWx < \\q\\x, < CohWx, M e qeX. 



Recall from Remark 2.1 that for parametric transport equations (3.35) is valid but (3.34) does not hold. 

At any rate, due to the denseness of X and Y, we can find sufficiently large but finite dimensional truth 
spaces Yj^ d Y, Xj^ C X, typically finite element spaces, that can provide a desired target accuray of the truth 
model. Since we are dealing here with problems for which standard tight a posteriori bounds are not available, 
we comment first on the truth certification. Note that this is particularly important for convection dominated 
convection diffusion equations when a complete resolution of very steep layers is prohibitively expensive even 
for the truth solution. We know that \\p{fJ,) — q\\x — 11/ ^ ^(j-IWy/^- order to be able to accurately evaluate 
the residual in the dual norm || • \\y^ one needs in any setting suitable assumptions on data oscillation, see e.g. 
[3 [311 [7]. One way to express this is to require that the projection of Ry^f into the test space Yf^ captures 
enough 

of i?yV- To this end, we make use of the following simple observation. 



Remark 3.8. Assume that (3.35) holds. Given W C X_\f and any S G (0, 1), there exists a finite dimensional 
test space V C Y such that 



M\\q-R^y;v\\^^<6\\q\\^^, 



q e Mx + W, fiEV, 



(3.36) 

which implies 

(3.37) (1 - _ B^Y' < WPY^yRY.if - B^.p)\\y, <\\f - B^Y', peW. 



In the following, we denote by V{W,6) all test spaces in Yj^ which satisfy the stability condition (3.36). 
Proof. Since A^x is compact there is a linear space Vm such that 



inf \\q-R-'Blv\\^ 



< d inf II g - 



q e Mx, M e ^■ 



It follows that the space Vm + W satisfies (|3.36|). Furthermore, since 
(3.38) 

vev 

and since inf„gy \\R^^ Bf^q — I'll = ||(/ — Py ,y)^^^ -S^fzlly,, the assertion follows. 



\nl\\q-R-^B^ 



\x. 



i^\\R-^B^q-v\ 



<5\\RfB^q 



□ 



We shall comment later how ( 3.36[ ) can be realized, see also [3 [32] for a more detailed discussion. Since 
rv,w{p, f) ■■= Py^^RyIU - Pfj-P) is given by 



(3.39) 



{RY,rv,w{p, f),z) = {f - B^p, z), z eYf^, 



the middle term in (3.37) is computable. 
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Remark 3.9. In what follows we shall always assume that for some fixed 5j\f < 1 and any given C X , the 



finite dimensional space Y_\f is contained in V{Xf^,S_\f) satisfying (3.36). Therefore, abbreviating the solution 
of ( |3.21[ ) for W = X_\f, V = Y^f, as PXm,Ym ^- PAfifJ^) G X_\f, u_\f{^) := UYj^^Xj^ifJ^), we immediately conclude 
that 



(3.40) 



\M~pM\\x..<i^-sT'^'\Mt^)\\ 



Remark 3.10. (i) For any desired target tolerance t, as soon as the computable quantity \\uj\f{p)\\Y^ drops 
below (1 — (5^)^^/^T we know that the truth solution has guaranteed accuracy < r which can be achieved by the 
refinement scheme in [7]. 



(ii) The above choice ofYj^f guarantees, by Propositions 3.6. 3.7 in particular, that 
h^,{q^v^ > i\ll-5l =: /3aA > 0, ^l e V, where C := ( ^' 



(3.41) 



inf sup 



lex^veYj^ \\v\\Yjq\\x,, 



j3j^ > 0, /i G where ^ 



and where (3 is the inf-sup constant from (3.31). Hence, j3j^ can, in principle be driven as close as one wishes 
to one or (3, depending on the choice of norm for X^. 

Note that the above statements do not contradict the possible case that the norms || • \\y^ or || • \\xi_,, eV, 
are not equivalent to a single reference norm. 

In the following, we shall often not distinguish for simplicity of exposition between truth and full spaces 
unless explicitly stated. In particular, whenever we speak of a computation in Yf^,X^ wc refer to the truth 
spaces endowed with the norms || • ||y^, || • ||x : respectively. 

Finally, the way how the bilinear forms depend on /i is important for practical feasibility. We assume that 
the dependence of the bilinear forms on /i is affine in the usual sense, i.e. 



(3.42) 



fc=l 



with parameter independent bilinear forms 6^, (•,•), fc = 1, . 



,mB, and smooth functions 0^ . 



4. Stabilization 



Suppose we are given a pair W = Xn C Xj^, Yn C Xv of finite dimensional spaces with bases <!>„ = {(f>j}^=i 
and ^I'n = {ipj }™f , respectively. Our convention will always be that the index n reflects the dimension of X„ 
while generally dimy„ — m{n) > n. While the purpose of Xn is to approximate Mx the role of Yn is, in view 
of Proposition |3.6[ to guarantee uniform inf-sup stability. More precisely, whenever Yn is (5-proximal for X„ 
( |3.22 1 for some (5 < 1, one has 

MyJWWx 



(4.1) 



inf sup 



> Vi - s^, iJ-(^v. 



Hence, a natural strategy is to choose a constant < C < 1, replace the right hand side of (4.1) by C,^/Y^^^ 



and enrich the space y„ until this relaxed inf-sup condition is valid. The closer one whishes Px„,Yn{l^) to be 
to the best X^-approximation pip), the closer C, should be chosen to one, see (3.231. In particular, any 

^ < 1 is in principle feasible. 

We shall formulate actually two variants of such a stabilization scheme which apply under slightly different 
assumptions. 

4.1. Inf-sup stabilization. The first natural idea which has already been used in is to enrich Yn 

by the supremizer for the infimizing parameter p.. More precisely, we first search for a parameter p, (z V and a 



function q £ Xn for which the inf-sup condition (3.27) is worst, i.e. 

bp.{q,v) . r ( . f bfj,{q,v) 
sup 17-|| |7:qi = mi mf sup -—r „ — 



(4.2) 
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If this worst case inf-sup constant does not exceed yet a desired uniform lower bound, F„ does not contain an 
effective supremizer for fl,q, yet. However, since the truth space satisfies the uniform inf-sup condition (3.41) 
there exists a good supremizer in the truth space which, on account of Remark [XT] is given by 

n-io - bp{q,v) 
v = RyBp^q^ argmax — , 

v€Yp \\v\\Yjq\\Xp 

and provides the enrichment 

(4.3) y„ -> span{K„, Ry^Bf^q}. 

This strategy can now be apphed recursively until we reach a satisfactory uniform inf-sup condition for the 
reduced spaces. 

Of course, three questions immediately arise: 

(i) Is the computation of fl and q feasible? 

(ii) Does this process terminate after finitely many steps? 

(iii) If so, is the number of necessary stabilization steps affordable? 

Assuming for the moment to have positive answers to (ii) and (iii), we first derive a suitable offline/online 
strategy for an efficient implementation of (4.3 1. First note that for given /I and q the new test function 
V := Ry. Bp^q can be computed by a standard Galcrkin scheme 

{"^^v^y- = bp,{q,v), V e Yp, 



so that it remains to solve the optimization problem (4.2) to find ft and q. To this end, we first rewrite the 
inf-sup condition in terms of the coefficient vectors with respect to the reduced bases. To describe this, we 
denote the corresponding Gramians, respectively cross- Gramians as 

((0»>j)xj^-^i, 

(5^(<?!)i,Vj))j.;"i- 

Practical feasibility relies on the following 

Assumption 4.1. In addition to B^ the Riesz maps Ry^ and Rx^ depend affinely on the parameter fi. 







■■= (*,*)y. 


(4.4) 













Remark 4.2. Under Assumption 4-1 the matrices in (4.4) can be computed online. Thus, by rewriting the 
left hand side of the inf-sup condition as 

mi sup -j--^ — 

(vTR^^vY^' (q^Rx,g)'/' 

we are left for each parameter /i with an optimization problem only of the size of the dimensions m{n),n ofYn 
and Xn, respectively. 

In order to find an infimizing q G M" we eliminate the discrete Riesz maps in the denominator by factoring 
them as 

(4.5) Ry, = L^^Ly^ Rx, - Lj.^Lx,. 

1 /2 

Here, one can think of a Cholesky factorization or of a spectral decomposition Ly^ = A-^ Qx^, where the 
columns of Qx^ form an eigenbasis and Ax^ is the diagonal matrix with the eigenvalues (in descending order) 
on the diagonal. Replacing v by Lh^v and q by Lx^q and defining := Ljj^BpL^^^ we find that 

u a\ ■ f v^B^q . v^D^q 

(4.D) mi sup , T-rz = ml sup 

qeR-^gK-i") (^;TRy^^)^ (g'TRx^q)^ ^gR™(") 

and hence, one easily verifies the following fact. 
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Remark 4.3. For any given the corresponding inf-sup constant is the smallest singular value of Z?^ and 
the optimal q is the corresponding right singular vector. Since the computational cost of the singular value 
decomposition is polynomial in the dimensions of the reduced bases, we can afford to compute all the inf-sup 
constants for a sufficiently large sample set S d V of parameters, yielding the optimal p.. 



The complete scheme is summarized in Algorithm [2] which we formulate for the general norms in (3.2) 



and the inf-sup constant /3jv from (3.41). 



Algorithm 2 Update to achieve inf-sup stability. 

1: function UPDATE-INF-SUP(y„, Xn) 

2: Choose < C < 1- 

3: Select a sufficiently large sample S C V. 

4: repeat 

5: for /i e 5 do 

6: Assemble the Gramians and cross- Gramians Ry^, , R-x^,, B^. 

7: Compute the Cholesky decompositions 

Ry — Lv Ly Rx = l'Z' Lx ■ 

8: Determine the smallest singular value (T(/i) and corresponding 

right singular vector of the matrix — Ly^ B^Lx^^ ■ 
9: end for 

10: Set p, = min{cr(/i) : € 5} 

11: Update Yn ^ span{Yn, Ry^Bf^qp,} with qp, = I]"=i(9/i)»</'i- 

12: until a{fl) > (I3j^ 
13: return y„ 
14: end function 



4.2. Stabilization based on (5-Proximality. We shall now formulate an alternative stabilizing scheme. It is 
related to greedy approximation and will shed some light on the above stabilization algorithms regarding the 



questions (ii), (iii). The idea is to enrich the space Yn to obtain stability based on the equivalent criterion (3.26) 
which can be rephrased as 

(4.7) inf WRy'Bpq - 011 Y^ < S\\Ry'B^q\\Y^, V q G X„, G T'. 
Defining 

(4.8) Xlip) := [q G X„ : H^H = pY.'^M'Zlk. = l} . 
this is equivalent to 

(4.9) sup sup mi \\Ry^Bpq-^Y,<5. 

We can again employ a greedy strategy to search for the parameter /i e "P and the element in Xlj^ := IJ^e-p ^nil^) 
for which the error is worst: 

(4.10) ifi.q)= argmax ini WRy^B^q - 

As long as the approximation error for p, and q exceeds some fixed 6 G (0, 1), we add the best approximation 
from the full truth space to the reduced basis: 

Yn spa.n{Yn, Ry^Bfiq}. 
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Since, as pointed out below (3.6), Ry^ = B^*R-^ B^^,we see that (Ry^y^ = Bf^R^^B^. Hence, in view of 
( |3.6[ ), we conclude that 

(4.11) argmax ( inf ||g - i? ) , 

which implies the following observation. 



Remark 4.4. // (3.35) holds so that all the spaces agree with a parameter independent reference space X, 
the output {p,,q) is the result of a greedy approximation to the compact set X\ := U,tg-p ^k(m)- Therefore, in 
principle, the scheme fits into the standard greedy theory in [51 [Tl [TT]. In fact, by ( |3.35 ), (3.42), and the fact 
that B* : — > (X^)' is an isometry, the set R^^B*Yn is a finite dimensional subspace of X . 



It remains to find a fast algorithm for the solution of the maximization problem (4.10) which will make use 
of the II • 11^ -norm (|3.6[) for X^. 



Lemma 4.5. Let q — Y^^=i1j't'j 
R-^ ;= , one has 



(4.12) 



q^$. Referring to the matrices B^,Ry^ from (4.4), and defining 
{I PY^.Y^RylB.qfy^ = (R^^ - 31^^)3,)^. 



l2 



Proof. By orthogonality of Py^.Yru &iid (3.6), we have 

W - PY^^Y^RvlB^^qfy^ - Ikll^^ - \\PY,^.Y,.RylB,q\\ 

By definition, we have ||(7||^ = q^Rj^ q. As for the second term, note that Py^^y^Ry^Bfj^q is the Galerkin 
solution of Ry^z = B^q. Since for any w G the coefficient vector z of Py^,y„Ry^w is given by 

Ry,z = (w,*) ((w,V'j))™i, 
we conclude that for w := q^B^^ one has 

z = Ry^ (*, B^*)q = Rv'S^q. 

Hence 

\\Py.,y.RY!:Bf.q\\Y^ = {B^q,Py^,y„Ry'^B^q)y^ = q^BjR^lB^q. 
which confirms the claim. □ 



Similarly, by the definition (3.6) of the X^-norm, we have \\q\\-^ = ||i?y^-B^(7|| so that 

X^ifi) = {qeX„: q^ q^*, Q^R-x^q = l} • 



It follows that the optimization problem (4.10) is equivalent to 
(4.13) 



(M,q) 



argmax 



qTR^ q 

where q is the coefficient vector of q. This problem can be solved analogously to the corresponding optimization 
problem (4.2 ), (4.2 ) of the inf-sup condition so that we obtain the alternative algorithm Update-(5 for updating 

The efficient practical execution of Algorithm Update-(5 requires assembling the matrices R-jj- in the typical 
offiinc/onlinc fashion. This is possible when instead of Assumption 4.1 the following holds. 

Assumption 4.6. The Riesz maps Ry^^, Rj^ and hence their inner products {■,-)y^, {'t')x depend affinely 
on the parameter fj, CzV. 
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Algorithm 3 Update to achieve (5-proxiniality. 

1: function Update-5(Y„, X„) 
Choose < (5 < 1. 

Select a sufficiently large sample S C V. 
repeat 



Assemble the Gramians and cross- Gramians Ry^, R-x ' -^m i^^^ (4.4)). 
Compute 



max 



q^(Rx -s^RYi^M)q 



{H, q) = argmax ^" -rp ^ 

Update r„ ^ span{r„,i?y_^Spqp} with 9^ = I]"=i(9/i)»'/'i- 

until 6.max < 5 

return y„ 
end function 



By (3.6), Assumption 4.6 is valid if the F^-norm can be chosen independent of /i, i.e., when (3.34) holds. 
Moreover Assumption 4.6 can also be satisfied for parameter dependent y^-norms as e.g., in view of (3.8), for 
the transport equation. 

Finally, it is important to note that the number of operations used by both algorithms Update-Inf-Sup 
and Update-5 (under Assumption 4.6) only depends on the size of the sample S and the dimensions n and 
m{n) of the reduced bases. Especially, it is independent of the dimension of the truth spaces which renders 
these algorithms feasible. 

Assumption 4.6 is clearly more restrictive, i.e. the use of Update-(5 is more constrained than Update-Inf- 
SuP which applies under the standard assumptions of afhne dependence and for any norm on A^. 

4.3. Interrelation between both Stabilization Schemes. We discuss next the interrelation between the 
schemes Update-Inf-Sup and Update-J. 



Proposition 4.7. Assume that we use the || • -norm for and the spectral decomposition 



in (4.5) for the scheme Update-INF-SUP, where is the diagonal matrix with the eigenvalues and the 
matrix of corresponding eigenvectors. Then the outputs of Update-inf-SUP and Update-(5 coincide. 



Proof Let M, 



BIR^^B^. Clearly, since is orthogonal, 



Amax(M) := max 



q^(Rx„ 



]_/2 ~r 1/2 

is the largest eigenvalue of the matrix I — A^ ' Q^M^Q^A^ ' so that 
(4.14) A.„ax(/i) - 1 - A,,in(A-i/2Q^M^Q>-i/2) 



On the other hand, using the A^-norm in Update-Inf-Sup, i.e., replacing Rx^ by R-jj- in (4.5), and using 
the spectral decomposition Rj^ = QjA^^Ay^Q^ for L^^^Lx^, the matrix in (4.6) takes the form = 
L^JB^Q^A^-*-/^. Clearly, the smallest singular value of is just A,„in(D^D^)-'-/^ and the corresponding 
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eigenvector agrees with the right singular vector of D^. Since 



we see that in this case the enrichments procued by both schemes agree, which confirms the claim. 
4.4. Termination of Stabilizing Greedy Loops. 



□ 



4.4.1. The general Case. Under the most general assumptions, neither insisting on (3.34 ) nor on (3.35 1 we resort 
to a very crude argument that ensures termination of the stabilizations loops Update-inf-SUP and Update-(5. 
Our findings can be summarized as follows. 

Proposition 4.8. Both schemes Update-Inf-Sup and Update-(5 always terminate after finitely many steps. 

Proof. We prove the assertion only for the scheme Update-Inf-Sup. The argument for Update-(5 is identical. 
To this end, let Yn+i and Xn+i be the spaces obtained by applying Update-Inf-Sup to the input spaces Yn 
and Xn. According to the update rule (4.3) used by Update-Inf-Sup, the enlarged space Yn+i is contained in 
the truth space Ya/- Thus, since Yj^ is finite dimensional, the statement of the proposition follows if each added 
function is linearly independent to the previous ones. To this end, assume the algorithm has already grown Yn 
to Yn and let v = Ry^Bp^q be the next function to be added (see (3.3|). Now, assume by contradiction that it 

is already contained in Y^- Since ■D is a supremizer this implies that 



sup 

veYn 



Recalling that p, and q are the worst possible choices according to (4.2), this violates the stopping criterion in 
Line 12 of Update-Inf-Sup. Thus, it follows that v is linearly independent from Yn showing finite termination 
of Update-Inf-Sup. □ 

The fact that, by the above argument, the number of stabilization steps may depend on the dimension of 
the truth space is certainly very pessimistic and not satisfactory from a practical point of view. In fact, much 
more can be said under some additional assumptions. 

4.5. Uniformly Equivalent Norms. Suppose now that all the spaces F^, agree as sets with Y, X, 
respectively (see ( 3.32[ )), and that (3.34), (3.35)) hold. Then we can replace || • ||y^ by an uniformly equivalent 
reference norm || • ||y . Since the Riesz map Ry is now independent of fj,, Assumption 4.6 holds and the stabilizing 
schemes Update-inf-SUP and Update-(5 are equivalent, see Proposition |4.7| Moreover, recall that, by Remark 
(3.2), the supremizer for g„ £ Xn in the inf-sup condition is given by Ry B^qn. The key observation is that 
because of the affine decomposition (3.42) of all these supremizers together generate a finite dimensional 
space. 



Remark 4.9. Given Xn C X, Xn = span{0j : j ~ 1, 



i}, let 



(4.15) Yn := {Ry'B^p : p e X„, fi e V} C span {Ry^B^^^ : j ^ 1, . . . ,n, fi e V} . 

Then Yn is a finite dimensional .space of dimension dimy„ < mgn, where mg is the number of terms in the 
affine expansion (3.42). Hence, one has 

(4.16) 



inf inf sup 



bf,iq,v) 



inf inf sup 



bp.{q,v) 



t^evqex^^^y^ \\q\\xJv\\Y t^ev qex„ ^^Yj^ MxJMy 



where fi^ > is the inf-sup constant from (3.41) in Remark 3.10 



Proof. If Bk is the operator corresponding to the bilinear form &fc(-, •) in the affine expansion (3.42) and Ry 
Ry is independent of e T', we conclude that 



Yn C span [Ry^BkC, 



J 



1, 



. ,n, fc = 1, . . . ,mb} 
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which proves the first part of the claim. Since all optimal test functions are contained in Yn the discrete inf-sup 
condition (4.16) follows immediately from the assumed inf-sup condition (3.41) of the full problem. □ 



The following simple observation is an immediate consequence of Remark 4.9 



Proposition 4.10. Assume that (3.34), hold. Then the update algorithm Update-Inf-Sup, and hence likewise 



Update-(5, increases the dimension of the test space in each step and terminates with a test space of dimension 
at most nniB- 

Proof. The proof is identical to the one of Theorem |4.8| by noting that all supremizers that are added during 
the algorithm are not only contained in the truth space Y_\f but in the much smaller space Yn C Ya/ which is of 
dimension msn. □ 

The above reasoning applies verbally to other saddle point problems like those appearing in parameter 
dependent Stokes systems or constrained optimization problems. The finite dimensionality of Y^ is also the 
basis of the a priori choice of stabilizers in [131 ttH [IS] to guarantee inf-sup stability although the connection 
with a greedy stabilization does not seem to be made there. 

The reason for nevertheless applying such a greedy stabilization is that a sufficient inf-sup stability might 
actually be achieved at an earlier stage so that in total fewer stabilizers suffice. 

4.6. A Greedy Perspective. As we shall see in later applications, in the context of Section |3.1| it will be 
important to treat also the case where only (3.35) holds but (3.34) is not valid. In this case the norms || • Hx^,, || • 

For instance, in 

IX, l^ eV, 



II are all equivalent and can be replaced by a parameter independent reference norm || • \\x . 
the case (|2.11|) one has || • ||L2(f2) = II • II ^ = W ■ \\x, ^ ^V, which will be further discussed in later numerical 



experiments. For the remainder of this section we assume that only (3.35) is valid. 

We have already seen that (9, /t) := argmin^g_Yi ( sup„gy^ w)/||f ||y^J agrees with the output of 



(4.10) and, on account of Remark 4.4 of a greedy approximation step to the set X^. Hence, the question of 



termination of the stabilization loop is equivalent to finding the smallest j for which 
(4.17) (T„^ 



max 

qex\4iev 



( ■mi \\q-Rl,^B;i,\\x)<5, 



where Y^ = Yn-i and Y^ is the enrichment of Y^ produced by the jth stabilization step. Here we assume that 
for the preceding pair (A„_i,y„_i) we have that K„_i = Y^i-\^ satisfies crn_i,£„_i < 5. We wish to see now 
how Y^ evolves from Yn-i- For convenience let Kf, Rx^B*^ 

A straightforward application of the currently available greedy concepts from [TJ [TT] is complicated by the 
fact that the sets become "less compact" when n grows and that the approximating subspaces R'j}B*^Yn 
depend on /x through the application of i3*. The following discussion is merely to shed some light on the 
expected behavior of a^n,j^ in particular, to identify some driving mechanisms, while we postpone a more 
detailed discussion to forthcoming work. 

Our first remarks concern the continuity of the mapping /i 1— >• K^. To this end, recall that the space 

Y 

individual norms 



Pl^gP Y^ is endowed with the norm 



from (3.33) which is here allowed to be stronger than the 



Y^. In view of (3.35) and (3.7), we have for any if! &Y 

WK^^nx < CoCmWK^^jW^^ - CoCmUWy, < CoCmMy- 

Thus, Kf^ G L{Y,X) which is equivalent to saying i?* G L(Y,X'). Now let B^. be the component of i3* 
corresponding to the fcth bilinear form 6fc(-, •) in (3.42) which, by assumption, are smooth. Obviously, one has 



(4.18) 



\{K,,-K^,)i:\\x<J2\ei{fi) 
fc=i 



< C max 

k — l,....ini 



\Ql{p^)-QW)\U\ 



Y- 



which shows that the mapping V — t- L{Y^ X), /i 1-^ is continuous in /z. By compactness of V , we can find for 
each e > a finite e-net comprised of Ne{V) centers n^j such that for each ip £ Y and any ^ £ V there exists 
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a J e {1, . 
(4.19) 



■ ,N,{V)} such that 



\\{K^-K,.Aij\\x<eM 



Y- 



In order to estimate HV'IIf we introduce the constant 
(4.20) C{n,Af) 



max 



\K-'q\\Y, 



which is finite because Kr, 



e IV C Y. However, a point of concern is that the quantity C{n,J\f) may depend 
on the truth space dimension, a point that wiU be taken up later again. In particular, we have for any given 
n and any ^ = K'p^^ipr, r < n, where (fir is any of the orthonormalized reduced basis functions generated by 
snapshots from M.x^ 



(4.21) 



\\{K-^-K^^^)K^^<j3r\\x < e\\K^'<j)r\\Y < eC(n,AA), fieV. 



For n = 1 the greedy stabilization would determine a sequence /ii,/, ^ = 1, , 



such that 



(4.22) 



sup 



E 

1=1 



< 6. 



X 



Thus, for e < 6 /2C{n,N), this means that £i < N^{V). In fact, as long as (4.22) does not hold no two fiij can 
fall into a single ball of the e-cover of V and as soon as every ball contains a jlu (4.21 ) says that (4.22) is valid. 

It is now easy to display for any given n a space Vn 'ZY which is (5-proxinial for X„. In fact, pick an e„-net 
for V where e„ := 6 / {C{n,N)^n) and let 

K := span{if^/^^0fc : / = 1, . . . , N,^{V), k = l,...,n}. 



Hence, for any ii £ V there exists a center such that, on account of (4.21), for every r 

110. - K^{K;^l jr)\\x = ||(if^,„,, - K^){K-Ur)\\x < 
Since the 4>k are X-orthonormal we obtain for any q e 



1, 



inf \\q- 



KM 



X 



< 



r=l 



ti 

\{q,cj,r)x\Ur - K^{K-l Jr)\\x < 5{j2\^'l^'t'-)^\ 

fc=l 



1/2 



< 5. 



Of course, it is not clear whether the greedy procedure for building the spaces Y^^ would actually produce a 
space of similar dimension dimT4 < nN^^{V) ~: Nn- If one did instead a separate greedy procedure for each 
subspace spanned by the argument for (4.22) would say that it terminates after at most 7V„ steps. Intuitively, 
one expects that the actual greedy algorithm terminates earlier since each individual 0. has to be resolved only 
with accuracy 6, not with accuracy e„ = 5 / {■n}/'^C{n,M)) as above. 

The perhaps simplest, although grossly pessimistic, way of rigorously bounding the number of greedy steps 
providing the spaces Y^ in the stabilization loop, is to use the above pidgeon hole principle and consider in 
addition to a ^-net for V with centers /^j Z = 1, . . . , Ng{'P), also a e-net for the compact set Xlj^ with centers 

I = 1, . . . , A^e(X^), where 5, e will be specified later. Now suppose that K^^q is the new snapshot added to 
Yji^^ to generate Y^. There exist, by construction, indices I G {1, . . . , Ng{V)} and i £ {1, . . . , Ni{X^^)} such 
that \\q-q^\\x < e and \\{Kf, - K^,^ ^)iP\\x < S^Wy- 



Then, for any (g,^) £ X^ x V with \\qi — q||x < e and ||(if^ 
to the new snapshot K^^q as above, we obtain 

W-K^ip\\x < \\q 



-f^^Ms.JV-'lix < Sll'iplW where q^,^J,s^l are related 



inf 



K.K^'qWx < 



q-q^\\x + \\q^-q\\x + WiKf. - Kf,)K^ ^q\\x 



(4.23) < 2{e + 6C{n,Af)). 

Now choose e and S such that 2(1+ 6C{n,J\f)) = S, e.g. by taking e 
(4.231 that a new snapshot Kj^ q can only satisfy inf^gyj-i \\q 

q\\x<e,\\iK^ 



{(9,m) : II* 



(5/4, 6 = 6/{4:C{n,U)). It follows from 
KptpWx > (5 if it falls into a cover element 
Kfig ,)4'\\x < '^llV'lly} that does not contain any previous snapshot 
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yet. This can happen at most Ng{'P)Ng{X^^) times which therefore bounds the number of possible greedy steps 
in the stabilization loop. 

As mentioned before, this bound is very pessimistic. In fact, since is isometrically isomorphic to a unit 
sphere in £3 the covering numbers Ng{X^) increase like (12/e)", see [Ml Chapter 13]. The numbers N^{V) 
instead depend only on the fixed dimension of the parameter set V and the smoothness of the parameter 
functions 6^(/x). 

One way to ameliorate the strong dependence of the 7Vg(X;^) on n is to relate the problem to a greedy 
approximation to a compact set that is independent of n. To this end, recall the solution set = Mx x {0}, 



see (3.181, which under the present assumptions is compact in X x {0}, independent of the truth spaces. As 
detailed later the spaces X„ are generated by a (weak) greedy algorithm. By compactness, the (weak) greedy 
errors 

(4.24) cr„ := cr„(A^x,^n) := maxdist(Alx,A:„)x ^ 0, 7W 00, 

tend to zero at a rate that is independent of the truth dimension. A repeated greedy approximation generates 
an AT-orthonormal system {(j>j}f^i C X. Let X° denote the closure of the span of {4>j}'jLi, i-e. 



3 en 



Let 



T:={qeX° : \qU<^}, sup ( y]{q,^j)^ 




Obviously, B := {q Cz J- : max{||(7||js:, < 1} is a compact subset of X and by construction 

(4.25) dist(g,A:„)x < cr„|(j|*, g e J", n e N. 

Moreover, the greedy errors for B are comparable to the greedy errors for Mx- In particular, 

(4.26) maxdist (6,X„)x< maxdist (TWjf, A:„)jf, n e N. 
Furthermore, 

Xi := {q € X„ : ||q|| . < 1, E V} C {q e X„ : \\q\\x < c^/} C T, 



since for q £ 

(4.27) |<Z|, -maxa-i( ^ (<Z,0,)2,) < <j;^'\\q\\x < c,, 



Therefore, recalling (4.17), we conclude that 



cr„j < max ( inf ||<7 - i?/ V'H^f )(cAfCr«) ^ 
(4.28) < max ( inf \\q - R^^ B;i,\\x)CMicM'Jn)-\ 



where cm, Cm are the constants from (3.6|. Hence, termination of the stabilization loop reduces to analyzing the 
the necessary number of steps needed to enrich = until maXqgxine./iSP (inf^gyj R^c < 

cmc^u/Cm- Clearly, N^{X}^ n S) < N^{B) where B is now a fixed compact set. We can now apply the same rea- 
soning as above with X^, e, 5 replaced by B, ICM<^nlCMt ^cmCu/Cm, respectively. This leads to the alternative 
bound Ng^^cM /(iCM){^)^SanCM /(CMC(n.N')){'P) for the maximal number of greedy steps. Note that in this case 
C{n,Af) can be replaced by 

C{B,Af):^ max \\K;\\\y. 

Since every q E X^ or q G S is a linear combination of snapshots B^^^f = p(/i;) and since the enrichments 
of the test spaces are of the form B^*Rxqj, qj £ X^^, they are linear combinations of elements of the 
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form B^*RxB^,}f. Since the operators B^*RxB^,} at least preserve the regularity of / the quantities ^q, 
q G {q ^ B), where now the inversion is understood in the infinite dimensional spaces, possess the required 
additional regularity in Y when / is sufficiently regular, see the discussion of the transport problem in Section 



Summary 4.11. We can now summarize the above findings as follows: 

(f) // the constants C{n,J\f) or C{B,J\f) are uniformly bounded independently of the choice of the truth 
spaces the stabilization loops terminate after a number of steps that is independent of the truth spaces. 
Their dependence on n can be bounded in terms of the metric entropy of X\ or the metric entropies of 
B and V , coupled in the latter case with the greedy errors cr^ for Aix- 

(2) The constants C{n,J\f), C(B,N) remain independent of the truth spaces when f is sufficiently regular. 

5. A Double Greedy Scheme 
We shall discuss now a greedy strategy for constructing reduced spaces X„,y„ for the saddle point problem 



(3.161 which is a weak formulation of (3.1). 



The basic outline of such a strategy looks as follows: 

• Stabilization: Given a pair F„,X„, enrich y„ until liY^,x,Sl^) — C/^7Vj M G where (3j\f is given by 
(3.41) and Pr^.x^ifJ^), /i S 7^, is the inf-sup constant (3.31) for the reduced spaces X„ and y„. 

• Approximation update: In view of the best approximation property (3.23), (3.24), we then improve the 
accuracy of the reduced spaces with the aid of a greedy step. 

That last greedy step, in turn, requires a tight residual based surrogate as detailed next. 

5.1. Tight Surrogates. Suppose now that the pair of spaces X„ C X, Yn C Y satisfy the (5-proximality 



condition (3.22 ) for some 6 G (0, 1) a nd abbre viate the corresponding solutions of (3.21 ) as PnifJ-) ■= Px„,y„ (m), 
u„(/i) := UY^,x„{^J■)■ By Propositions 



3.6 



3.7 



the definition (3.6) of the X^-norm says then that 

IIp(m) -Pn(Ai)lljf„ = \\f - B^p„{fi)\\Y', tieV, 



i.e., the residual based surrogate 
(5.1) 



Xn X y„ 



1/ - s^p„(m)I|y; 

is in this case almost ideal. In fact, combining with (3.23) yields 

1 



(5.2) 



qy < R{ii,X„ X y„) < 



1- S qex. 



inf IIp(^) 



Hence, ( 1.8 ) holds with cn = 1 — S, Cn = 1. 



5.1.1. Reduction to Truth-Riesz Maps. Of course, the dual norm 
computed exactly. Instead, defining 



and hence R{n,Xn x F„) cannot be 



(5.3) II • 

we consider the following candidate 
(5.4) 



IP^ 



veY' \\v\\y. 



Rn{tj) := 11/ - S^P„(m)||y^, 



where we continue to assume that Yj^f e V{X^,5^)^ i.e., the truth spaces Xj\f,Yj\f comply with Remark 3.9 



Then, by (|3.39|) and (|3.40|), we conclude that 



My^-PrMx, < (l-^lA)-^/'||PY„Y^(i?y;(/-S^P„(M))l|Y. -(l-'5l/-)"'/'i?«(M) 
(5.5) < (1 - 5X-)-V2||/ _ B,pr.{^,)\\y, < (1 - S'^r'/'{1 - 5)-l|b(/i) - Px,,xAl')\\ X, 

This immediately implies the following fact. 
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Proposition 5.1. Under the above assumptions on the truth spaces the surrogate RnifJ-) given by (5.4) is tight 
with condition 



(5.6) 



K{Rn) < 



1 



(l-<5X-)i/2(i-<5)' 



which, in principle, can be driven as close to one as one wishes, at a computational expense caused by a 
correspondingly large truth space Yj^ and a possibly larger number of stabilization steps. 



The equivalence 
(5.7) 



\\f-B,M 



Y'r 



\f - B^qWv', qeXj^, 



which is nothing but a reformulation of (3.37) for W — X^, Yv G V{X^, Sj\f), says that the || • ||y/ -norm yields 



still a meaningful error estimate even in case the truth spaces are not rich enough to resolve all features of the 
infinite dimensional exact solution which will be seen below in the experiments. 
The above findings can be summarized as follows. 



Proposition 5.2. If both, (3.34) holds, then Rn{fi) defined by (5.4) is feasible. 

Proof. Under the given assumptions the norms || • can be replaced by a uniformly equivalent equivalent 
reference norm || • ||y so that the Riesz map Ry is independent oi fi €V. Hence, RnifJ-) can, in view of (3.39), 
be efficiently evaluated by a standard offline/online decomposition, see e.g. [27] . 



□ 



5.1.2. Rerative Tightening. Recall that in the pure transport problem (3.34) does not hold, see Remark 2.1 



Hence, the surrogate Rn{^J) from (5.4) is no longer feasible in the strict sense. Instead a feasible variant would 
be 



(5.8) 



where the dual norm is now induced by the reduced space Yn instead of the truth space Yj^. While the 5- 
proximality of y„ for Xn (see (3.22)) does ensure the equivalence HS^qUy^; ^ ||B^(7||y^, q e X„, (with constants 
close to one, depending on 5) the analog is not clear for |j/ — S^q||y^ since generally / ^ B^{Xn). 

However, Remark 3.8 immediately tells us at least a criterion for the validity of the desired residual equiva- 
lence, namely with the aid of a somewhat strengthened (5-proximality. 

Remark 5.3. Assume that for some 6 G (0, 1) one chooses Yn G V(X„,5) so that 



(5.9) 

Then 
(5.10) 

and we 



inf \\p-R-^'b;v\\x,<S\\p\\x,, ypeMx+X„ 



veY„ 



il-py/^f-B.qWr;^ < Wf-B^qWr;^ < ||/-B^g||y; 
have n{R'n) < (1 - 52)-i/2(i _ 



Note that we could replace Mx in (5.9) by its truth approximation Mx,Ar since, in view of (5.7), it suffices 
to establish ||/ — B^gjly^ ~ 11/ — -B^^llr^- But the main practical issue remains how to find Yn satisfying (5.9) 
at affordable cost. 

To this end, we propose a systematic way of successively substantiating tightness of error estimators at the 
expense of an additional computational effort in the offline phase. We refer to this process as iterative tightening. 



The idea is that once a reduced space provides sufficiently accurate approximations to A^X: condition (5.9) 



becomes easier to fulfill. To make use of this observation, assume we have a second pair of reduced spaces 
X C Xj^ and Y C Yj^. We now describe how such spaces can give rise to tight surrogates and later discuss 
their construction. 
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Lemma 5.4. Assume that the pair Xn + X and Y satisfies the (standard) 5-proximality condition (3.22) and 
that the approximation of p{^) from X„ + X is superior to the approximation from Xn alone, i.e., one has for 
some < ^ < 1 

(5.11) ||p(m)-p"(a^)IIx, <eib(M)-p"(M)IU,, f^eV, 

where p"{fJ.) and p"{lJ.) are the respective best approximations to p{fi) from Xn and Xn + X . Then 

(5.12) M\\p~R^Ib;v\\x,<~S\\p\\x,, ypeMx+Xn,fieV, 



where 6 :— {l + 5)£,+6. Hence, for £^,6 sufficiently small, the surrogate R'nip) from (5.8 1 is tight with a condition 
given by (5.6) with Sj\f replaced by S. 



Note that we use the space X„ + X as opposed to X alone because for the latter space the condition (5.11) 
would imply that Xn C X if both spaces are constructed from snapshots, which would be too restrictive for the 
application below. 



Proof. Then, for each deviation p(/i) — p, p € Xn, we obtain 

m^\\p{^,)-p^Rl,]B;v\\x, < 

(5.13) < \M-p-ifi)\\x,^+S\\p-ifi) 



b(M) -p"(m)I1x, + inf ||p"(m) ~p-RxIb;v\\x, 



< {{l + S)^ + S)\\p{f,)-p\ 



Thus for (1 + (5)^ + S sufficiently smal l the extended (5-proximality condition (5.9) is satisfied for the trial space 
Xn and test space Y. Thus, Remark 5.3 applies which says that the surrogate -R^(/i) from (5.8) is tight with 
the claimed condition. □ 

We shall describe ways of constricting the spaces X and Y later in Section [g] 

5.2. Approximation Update. Either scheme Update-Inf-Sup or Update-J outputs a pair X„,y„ that 
is uniformly inf-sup stable, i.e., the corresponding inf-sup constant is uniformly bounded away from zero 
/3F„,x„(Ai) > C/3aA, H & V. By Proposition 



5.1 



the surrogate i?„(/i), defined by (5.4), is tight with a con- 
dition controlled by the (5-proximality parameters. The feasibility of this surrogate depends on the way how 
the spaces depend on the parameter fi d V, see Proposition 5.2 In applications, an infeasible surrogate is 
replaced by R'nifJ-) from (5.8) combined with iterative tightening. 



This suggest the following oMier greedy step update- approximation, defined in Algorithm 5.15 which aims 
at improving on the accuracy of the reduced model. 

5.3. Putting Things together. The overall double-greedy method (see Algorithm DG-1 below) for computing 
reduced spaces for (7.1) consists now in combining the inner greedy stabilization loop with the outer greedy 
approximation step for the saddle point formulation (3.16). 



To analyze of algorithm DG-1 recall the solution manifold A4 = Mx x {0} from (3.18). Since the inner 



stabilization loops ensure, by Proposition 3.6 and (5.2), tightness of the surrogates, we can invoke Theorem 1.3 
The above findings can now be summarized as follows. 



Theorem 5.5. Assume that (3.35) holds. Let pn{p) := ?'x„,y„(A'); u„(/i) := uy„.x„(m) denote the solution 
components of ( 3.21[) fo r W ~ Xn, V = Yn, were [y„, A„] are the reduced spaces produced by algorithm DG-1 
using the surrogate (|5.4|). Let 



<Jn{Mx) ■■= sup ||p(/i) -Pn{fJ^)\\x > 

pep 



dni-Mx) ■= inf ( maxdistfAl^, Z„)x 

diin(Z„)=ri V 



(a) Then, if dn{Mx) = 0{n "), for some a> or if dn{Mx) — 0{e 
(5.16) a„(A^x) - 0(n-"), a„(A^) = 0(6^™°), 



), for some c,a>0, one has 
n oo, 
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Algorithm 4 Update to improve the approximation quality. 

Require: Finite dimensional spaces Yn C Y, X„ C X that satisfy the inf-sup condition 

(5.14) inf inf sup ^^f%^ > C/^aA 

i^eVgex^^^Y^ \\q\\x, 

for some < C < 1- 
1: function Update-approximation(F„, X„) 
2: Select a sufficiently large sample S C V. 
3: compute 

jj. :~ argmax Rn{n) 

4: Compute the solution [u,p] G Y/i x Xfi of 

iu,v)Yi,+bf,{p,v) = veYf,, 
bf,{q,u) = {g,q), q E Xf,. 

5: Set 

span{X„,p} Xn 
If (3.35) holds orthonormalize {(fii, . . . ,0„-i,p„} in X. 
return Yn, X„ 
end function 



Algorithm 5 Double greedy scheme 

1: function DG-1 

2: Initialize Yi = {0} and Xi — span{u(/xi)} for an arbitrary /xi e V. 

3: Fi ^ Update-inf-sup(Fi, Xi). 

4: while max^gp Rn{lA > t do 

5: X„ Update- Approximation(F„, X„). 

6: F„ -S— UPDATE-INF-SUP(y„, 

7: end while 
8: return Ym Xn 

9: end function 



respectively, with constants depending on the parameters S_\f,5,(^ in Update-(5 or Update-Inf-Sup, ana 

( |5.16 ) remains valid for an{M.x) replaced by 

P7i{p)\\x + - Un{fl)\\Y^}. 



the constants in (3.35), (3.7) 
(5.17) 



Moreover, 

CTniMx, 



sup {\\p{p) 



(b) Assume that both (3.34| and (3.35) hold. Then, the assertion (a) holds where in addition dim (1^ x Xn) < 
{l + mB)n, n € N. All bounds remain valid up to the tolerance tol* when all computations are carried out within 
this accuracy. Moreover, the surrogate (5.4) in algorithm DG-1 is feasible. 



Proof. The output [y„, X.^] of step 6 in DG-1 is uniformly inf-sup stable so that the surrogate (5.4) used in step 
5 is uniformly tight, with a condition depending on the stabilization thresholds 6, 6_\f. Concerning an{Mx) we 
use (3.24|. By ( 3.35[ ) the surrogates remain uniformly tight for the reference norm || • Hence, Theorem 1.3 



applies. The rest of the assertion follows from Proposition |4. 10 and Proposition 5^ 



□ 



In general, under the assumption (a), the well conditioned surrogate (5.4) is not feasible. Employing the 



feasible surrogate (5.8) instead, requires, in order to gurarantee rate optimality, an additional iterative tightening 
as described in Section 5.1.2 and later in connection with numerical experiments. Note also that under the 
assumption (a) dimy„ could be significantly larger than n, see the discussion in Section 4.6 In the case of 
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uniformly equivalent norms, i.e., when both conditions (3.34) and (3.35) hold, the dimension of the stabilizing 



spaces Yn remains proportional to the dimension of the reduced primal space. 

We conclude this section with a remark on the online evaluation of Pn(/i). Recall that the corresponding 
component € Yn is only an auxiliary variable tending to zero. 



Remark 5.6. Assume that (3.34) and (3.35) hold, i.e. the spaces and can be choose parameter indepen- 



dent. Instead of solving for a given /i the saddle point problem (3.21 ) for W = Xn, V = Yn, whose dimension 
is n + m{n), one can compute in the offline phase the test basis functions 4'k.j, j — 1, . . . n, k = I, . . . , ms 



(5.18) 



V e Yn, i = 1,. 



. ,n, 



where bk are the components of the affine decomposition ( |3.42 I . Then, defining 

(5.19) ^,"(/i):-^e^(M)V'fe.,, 

fc=i 

PXn,Yn{p) of the saddle point problem 



on account of Proposition 3.5: for each fi £ V , the solution Pn{l^) 
( |3.21 1, also solves the Petrov-Galerkin problem 

(5.20) ^A.(Pn(M),V';(M)) = (/^,^;(A')>, J = 

Hence the online complexity is indeed determined by the size n of the trial basis. 

6. Application to the Model Problems 



6.1. Singularly Perturbed Convection-Diffusion Problems. We refer to the setting in Section [2.3.11 and 
consider the convection-diffusion problem (2.5) for large Peclet numbers. 



To this end, we shall briefly discuss two scenarios concerning the truth spaces, namely (a) boundary layers 
are to be resolved completely by the truth spaces, and (b) due to a possibly very small diffusion, even the truth 
spaces are not required to resolve boundary layers. 

In case (a) solutions in the truth spaces could be obtained by simple standard Galerkin discretizations and a 
modified variational formulation according to (3.6 ) is only needed for the computation of reduced basis functions 
which then also resolve boundary layers well. 



In this example we prescribe the space Yu and adjust according to (3.6). We first decompose i?u into its 



symmetric and skew-symmetric parts: 

1 



B^u, v) 



{B^v,u) 



and define 



(6.1) 



-div6(/i) 



1/2 



L2(n) 



see 1301 E] for details. \\ ■ \\y^ is then equivalent to the standard _ff^(fi)-norm with constants depending on 
the diffusion e. This works perfectly when the discretization (adaptive or not) resolves the boundary layers. 



However, when layers are not resolved, although stable, the scheme (3.21 ) would give rise to unpleasant numerical 



artifacts, due to the nature of the involved norms, see the detailed discussion in |B1 15^. 

Therefore, we briefiy recall next an alternative variational formulation of (2.6) avoiding the numerical arti- 



facts, regardless of choosing sufficiently large truth spaces that fully resolve boundary layers or not. In essence, 
in case the finite element truth space does not resolve boundary layers this scheme behaves like a solver of 
the corresponding transport problem for e = which is, however, ill-posed when insisting on zero boundary 
conditions on all of We resort to a remedy proposed in [HI [35]. We retain the construction of the norms 
but modify the outflow boundary condition. Instead of building them into the trial space, we impose them only 
weakly. To this end, let 

r+(/i) := {xedil: n{x) ■ 6(/i, x) > 0} 



28 



WOLFGANG DAHMEN, CHRISTIAN PLESKEN, AND GERRIT WELPER 



be the outflow boundary where n{x) is the outward unit normal at x. Now, we take 
(6.2) X^:^{qeH\n): q\r_(,,)^0} 

as a set with norm defined below. Here and in the following, restrictions to the boundary are implicitly 
considered in a trace sense. Thus, zero boundary conditions are only built into on part of the boundary. To 



find a weak form of the boundary conditions at the outflow boundary recall from (3.11) the connection of (3.6) 
with the optimization problem 

(6.3) \\f-B,,p 



2 



where p belongs now to the larger space X^. So far we have not changed which is still i?o(i^) endowed 



with the norm (6.1). Due to the missing outflow boundary conditions, has a nontrivial kernel so that the 
optimization problem is not uniquely solvable. One simple remedy is to add the outflow boundary condition as 
a penalty term: 



(6.4) 



where 



\\f-B^p\\l-,+oj\\p\\j,^^^) 



1 /2 

is a norm for Hqq (r+(/x)) and w > 0, see |Q. Practically, this weak enforcement of the outflow 
boundary condition applied to a subspace W C Xf^ has the following effect: typically boundary layers are found 
at the outflow boundary which are too narrow to be resolved at affordable cost. If u is chosen small, then 
the enforcement of the outflow boundary condition has little weight so that it is almost ignored which, in turn, 
removes layer artifacts. If, however, W is sufficiently rich so as to resolve layers, inig^w \\f ~ Bf_tq\\Y/^ becomes so 
small that the boundary penalty becomes important and the boundary conditions are approximately satisfied, 
see Figure [la[ The rationale is that as long as the layer is not resolved the error with respect to conventional 
norms (including the SUPG-norm) is mostly concentrated in the layer region, which therefore stays, roughly 
speaking, as large as not realizing the boundary conditions at the outflow boundary at all. Putting a small 
weight on this error contribution actually increases accuracy away from the outflow boundary, see [B]. Putting 
it in a slightly different way, by allowing more freedom in the outflow boundary layer, the X^-norm is changed 
in such a manner that the error in the boundary layer has very small weight. This in turn allows one to better 
control the error away form the layer. 

To apply the theory of Section |3.H we define the test space 



(6.5) 



\\[v,g]\\l := 11^111-+ 



and the operator 

B^p := [B^p,p\r^(f,)]. 
According to the definition of the graph- norm ( |3.6[ ), this yields the norm 

:= \\B^.p\\h = 



Ml 



\B^p\\y' 



for the trial space X,, 



Note that in the optimization problem (6.4) the norm of the boundary penalty is not a dual norm. This 



allows us to replace the system (3.21), which in or case is a 3 x 3 block system, by the simpler system 

- {B^pw.v) 

- f-i{pw,q)HtM 



(6.6) 



(i?Y- My, 



{B;uy,q) 



0, 



veV, 
qeW, 

in all practical computations. This system is derived by the same reasoning as in Section [XT] applied to the first 
|y; of the optimization problem (6.4) only. In |Bll32j. it is shown how to transform this system 



term\\f-Bf,pw\ 

to the equivalent saddle point problem (3.21), so that the theory of the present paper still applies to (6.6 1 



In summary, we have found a stable variational formulation of the convection-diffusion problem (2.6) that 

However, note that the spaces and may differ even as sets 



fits into the general framework of Section 3.1 



for different fj. see (6. 51, (6.2). Specifically, the dependence of on lies only in the boundary conditions. 
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However, for a polyhedral domain 17 one can find a finite cover {Vi : I = 1, . . . ,P} of V so that the outflow 
boundary portions r-|_(^) = T+j stay the same for ^ eVi- Hence, the spaces F^, all agree as sets for fi £Vi. 
Clearly, the solution manifold M (see (3.18)) is a finite union of solution manifolds corresponding to the 

subsets Vi- Since each A4{1) is compact so is the finite union A4. Note that for fi G Pi the Riesz map Ry^ is 
independent of fi. Therefore, we can apply Theorem |5.5| to each component Vi leading to the following result. 



Corollary 6.1. The scheme DG-1 based on (6.6) is rate-optimal. 



In this case the online evaluations can be based on Remark 5.6 Some first numerical experiments are 
presented in the following section. 



6.2. Numerical Experiments for Convection-DiflFusion Problems. 

problem 



We consider the convection-diffusion 



(6.7) 



eAp 



f cos /i 
I sin ji 



Vp + p=l, in = (0, 1)^ p = 0, on dn. 



In all test cases we use the variational fomulation based on (6.4) regardless of the choice of the truth spaces. 
First, we treat scenario (a), i.e., with e = 2^'^ which is already convection dominated. However, we use a truth 
space that completely resolves the layers. Specifically, for X_\f and Yj^ we choose bilinear finite elements which 
are continous on a rectangular uniform grid of meshsize 2~^ and 2~^", respectively. For all computations we 
used an equidistant sample set S C V = [0.2, tt — 0.2] of cardinality 500. Using finite element a-posteriori error 



estimators from [B], based on (5.1 ), the respective truth-accuracy is bounded by 0.00569966. We note that these 
a-posteriori bounds represent the truth residual and hence the energy error only within some fixed constants. 
This is in contrast to the surrogate bounds for the reduced spaces which are much tighter. The number of 
adaptively generated basis functions for the reduced test space together with the corresponding constant of 
the (5-proximality, as well as the maximal surrogate are given in Table 1 a) Figure la shows a reduced basis 



solution for the angle /i = 0.885115. The parameter dependent direction of the first order term is visualized by 
a plane. 



(a) e = 2 ^, maximal a-posteriori error 0.00569966 (b) e = 2 ''' , maximal a-posteriori error 0.0197011 



dimension 




maximal 


surr / 




dimension 




maximal 


surr / 


trial 


test 


6 


surrogate 


a-post 




trial 


test 


6 


surrogate 


a-post 


2 


3 


2.51e-01 


7.04e-02 


1.24e-H01 




2 


5 


8.92e-03 


1.25e-01 


6.37e+00 


3 


6 


3.74e-01 


3.08e-02 


5.40e-H00 




3 


8 


1.22e-01 


9.65e-02 


4.90e-f00 


4 


7 


3.74e-01 


7.43e-03 


1.30e+00 




4 


11 


1.13e-02 


3.21e-02 


1.63e-h00 


5 


10 


3.51e-01 


5.81e-03 


1.02e-t-00 




5 


14 


1.27e-02 


2.61e-02 


1.32e-h00 


6 


13 


1.86e-01 


5.70e-03 


l.OOe-t-00 




6 


17 


5.55e-03 


2.21e-02 


1.12e-l-00 












7 


20 


4.82e-03 


1.97e-02 


l.OOe+00 








(c) e = 2 maximal a-posteriori error 0.001055 






dimension 




maximal 


surr / 


dimension 




maximal 


surr / 


trial 


test 


S 


surrogate 


a-post 


trial 


test 


S 


surrogate 


a-post 


2 


5 


1.35e-03 


2.11C-01 


2.00C+02 


12 


33 


3.47e-04 


1.60e-02 


1.52C+01 


4 


9 


1.09e-02 


7.58e-02 


7.19e-(-01 


14 


39 


l.lOe-04 


8.46e-03 


8.02e+00 


6 


15 


1.61e-03 


5.02e-02 


4.76e-|-01 


16 


45 


9.39e-05 


7.87e-03 


7.46e+00 


8 


21 


7.99e-04 


2.39e-02 


2.26e-|-01 


18 


51 


6.11e-05 


7.69e-03 


7.29e-f00 


10 


27 


3.55e-04 


2.10e-02 


2.00e-|-01 


20 


57 


5.28e-05 


6.35e-03 


6.02e-F00 



Table 1. Numerical results for the convection-diffusion problem (6.7). 
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(a) 



{h)e: 




(c) e 



Figure 1. Solutions of the convection-diffusion problem ( |6.7[ ). (a) RB dimension n ~ 6, 
m(n) = 13, angle fj. = 0.885115, (b) RB dimension n = 7, m(n) = 20, angle ^ = 0.257484, (c) 
RB dimension n = 20, m{n) — 57, angle /z — 0.587137 



-26 




0.1 
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2 4 6 5 10 15 20 
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Figure 2. Surrogates of the reduced basis approximation for the convection-diffusion problem (6.7). 



The intermediate case e = 2 shown in Figure lb demonstrates how the formulation handles a not fully 
resolved boundary layer which is not far-off being resolved either. 



The other example, scenario (b), refers to the same problem (6.7) again, however, with a very small viscosity 
= 2-26 Hence, this case is even more strongly convection dominated and poses difhculties for resolving 



the boundary layers even for the truth space itself. Since the boundary layers are not resolved (6.4) does not 
strictly enforce strong boundary conditions at the outflow boundary even in the truth space. Accordingly, the 
approximate solutions from the reduced space do not satisfy the boundary conditions in a strict sense either. 
In fact, we choose the same truth spaces as in the preceding experiment. Thus we have to employ the norms 

The numerical results are summarized in Table 1 [c) and a 



(6.5) based on the variational formulation (6.6) 



corresponding reduced basis solution is displayed in Figure [Tc| 

Figure [2] displays a surrogate plot for the values 2~^, and 2-^6 of e. One observes that the error of 
the reduced basis approximation decays rapidly already for small reduced bases. In fact, since infpgx^ 11/ ^ 
B^pWy^ < infpgx„ \\f ^ Bfj_p\\Y^, the error of the truth approximation is always a lower bound for the error of the 
reduced basis approximation. This contrasts standard reduced basis methods where one, in our terminology, 
chooses Yf^ = Xj^. Instead, we assume a larger space Yj^ e V{Xj^,5j^) which, according to Remark 3.8 
implies that the surrogate ( |5.4[ ) is equivalent to the true error with respect to the infinite dimensional solution, 
regardless of whether the truth space resolves all solution features like boundary layers or not. Comparing the 
surrogate plots with the Tables [l| one sees that the reduced basis errors stagnate roughly at the error level of 
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the truth solution. Due to the very small J-proximality thresholds, the surrogates reflect the true reduced errors 
very accurately, see 



, Proposition 5.1 



Note that in all cases the inner stabilization loop produces at most ttib — 3 addidional test basis functions 
for the test space, see Proposition |4.10| 



6.3. Transport problems - the Worst Scenario. We address now the transport equation (2.8) in Section 
2.3.2 Aside from its essential appearance in more general kinetic models and Boltzmann type equations, it 



can be viewed as a "limit" of convection-diffusion problems. The particular interest lies in the complete lack 
of viscosity as a "classical" stabilizing ingredient, see e.g. |23]. Moreover, as we shall see, the conditions (3.35) 
and (3.34) do not hold simultaneously, not even for suitable subsets of V. Moreover, the parameter dependence 
will be seen to be significantly less smooth. 



We have already proposed a variational formulation (2.9) along with the spaces Y^,X^ in (2.10) endowed 
with the norms (2.11). With these definitions, the operator i?^ : — )■ Y"^ is an isomorphism with condition 
number 1, i.e. || ■ = !| ■ 1| jjf = !l ' 1Iy,;i see [7 . Notice that in this case the Riesz map i?^^ are given by 

Ry^ = -S^-B*, i.e., {v,w)y^ = {B*^v,B*^w}, 



\L2{n), 



so that (3.42) and (3.351 are valid. Finally, the Riesz maps Ry^ and Rx^ depend affinely on the parameter so 
that the double greedy scheme can be applied. 



However, since the y^-norm is not independent of /i, we cannot evaluate the surrogate given by (5.4) in the 
usual way. As a remedy, we use the surrogate from (5.8), i.e., we approximate this inverse Riesz map by 

projecting on the reduced basis space Yn instead of the truth space Y_\f. To ensure that this surrogate is also 



tight we take up the criterion in Remark |5.3[ Specifically, we wish to apply Lemma 5.4 and try to construct 
suitable pairs X, Y as follows. 

We run the double greedy scheme (possibly) several times which yields the sequences of reduced spaces 
XI, X2, ■ ■ ■ and Y]*, F2 J • ■ • J i = 0, 1, 2, .. in the ith run of the full double-greedy algorithm. Now, say we stop 
the first run at index N and define X := Xpj. For the second run, we use the same inital spaces as for the 
first run, however, the calls of UPDATE-lNF-Sup(y„\ AT^) are replaced by UPDATE-lNF-Sup(y„\ A -f A^), so 
that (S-proximality is guaranteed for the larger space A -I- Xlj^. Then, with the n-dependent choice Y = Y^, 
the estimate (5.13) implies that for the second run the surrogates are tight as long as the condition (5.111 is 
satisfied. 

Of course, neither can this latter condition be rigorously checked since we cannot rely on the surrogates, nor 
have we specified the terminating index N = Nq. We briefly sketch now several options of iteratively tightening 
the surrogates R'„{^). One could stop the first run i — a,t the smallest Nq for which R'^^{ijl) / tjs/ < a for some 
a ^ 1, where T_\f is the truth error tolerance. The second run i — 1 with A = A^ = X'^^ will stop at step Ni. 
In general, the ith run with A* — A*^^ -I- A]y^^^ stops at Ni. One expects that A^^+i > Ni since the surrogates, 
being lower bounds for the true residuals, become tighter as long as A* grows. A practical stopping criterion 
would be, for instance, that Ni+i < N„ or " Bf,pN,ifi)\\Yji^ ~ -Rwi+i (m)/II/ - -B^m+i (Ai)llYX^• 

An alternative strategy is to apply the double greedy scheme to the defect problem 

and form A as the sum of Ajy^ and the largest reduced space for the defect problem. Since the relative accuracy 
to be achieved for the defect problem only needs to meet the constant ^ in (5.11 1 one expects that a few steps 
suffice. Since A now contains "complementary" information A^ is enlarged more effectively than in the first 
method. 

The upshot of these comments is that investing additional computational offline effort is guaranteed to tighten 
the surrogates and thereby improves the choice of the reduced spaces. This is in contrast to greedy strategies 
based on surrogates that are not based on well-conditioned variational formulations and therefore most likely 
fail to detect the most effective snapshots. These issues will be addressed in forthcoming work 

Since the basis function (jjj € A„ can now be orthonormalized in L2(^) and || • Hj^ = 
applies and yields the following result. 



\L2{n} 



, Theorem 



5.5 
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Corollary 6.2. If R'^[^) from (5.8) is based on iterative tightening with satisfying (5.111 for sufficiently 

small S^, then the scheme DG-1 using Update-inf-SUP is rate- optimal for Mx ■ 



We could also reverse the roles of the spaces Y^^Xu, choosmg L2{^) as the test space, see [TH 



In this 



case the trial spaces would essentially depend on the parameter /i so that the understanding of the solution 
set is less clear. On the other hand, this choice would correspond to the limit of the formulation (6.6 1 for 
vanishing viscosity. 



Since (3.34) does not hold we cannot apply Proposition 4.10 to predict a strict a priori bound on the number 



of stabilization steps in Update-;? or Update-Inf-Sup. Adhering to the notation in Section |4r6| we have here 
Kfj^ = i?*, see (2.111. Since in the present case Rx is the identity, as pointed out there, the enrichments of 
the test spaces are are linear combinations of elements of the form BJ^*B^,}f where , are different 
most of the time, due to the greedy selection. As a consequence, when / £ L2(ri), this means that indeed 
B~*B~,}f S H^{n) = Y. Of course, the H^-norm may deteriorate when fi' , ij" get closer, which however 
may be offset to some extent by the expectation that these snapshots are most relevant for the stabilization 
of solutions with nearby parameters. It is also clear that higher regularity of / would indeed ensure sufficient 
regularity of the q S (or q £ B), independently of fi',fJ." and hence allows one to control the constants 
C{n,N) (or C{B,N)). This effect is reflected to some extent by the experiments below. 

6.4. Numerical Experiments for transport problems. We consider the analog of the convection- diffusion 



problem (6.7) with zero diffusion e = and corresponding boundary conditions, i.e. 



(6.8) 



/ cos /i 
\^sin/i 



Vp + p = 1, in r2 = (0, l)^ p = 0, on r_ 



We employ a truth trial space with mesh size 2~®, using discontinuous piecewise bilinear finite elements with 
proper boundary conditions. To ensure stable truth discretizations, the test truth space is comprised of globally 
continuous piecewise bilinear finite elements, therefore being contained in F = fl/jsp ^ finer mesh with 

mesh size to ensure 5-proximality. Recall that the spaces now differ even as sets. The results are shown 
in Table [2] and a reduced basis solution for the angle ^ = 0.244579 is given in Figure [3j Specifically, in addition 
to the dimensions of the trial and test spaces in columns (1) "trial", (2) "test", it records the values of the 
surrogates in column (4) "surr" , the error between the reduced basis solution and the best L2-approximation 
of the exact solution in the truth space in column (6) "rb L2", the error between the reduced basis solution 
and the truth solution in column (5) "rb truth", and finally in column (7) "surr/err" the ratio between the 
computed surrogate and the error in "rb L2" . All values reflect the worst case over the parameter range. 



As pointed out above, unlike the convection-diffusion problem, the surrogate (5.8) for the transport problem 



is not necessarily well-conditioned from the start. Therefore, Table [2] contains one column which shows the ratio 
of the surrogate compared to the true error of the reduced basis approximation, maximized over a sample of the 
angles with the largest values of the surrogate. Although this is at this point not founded rigorously, we see that 
this ratio stays uniformly bounded with respect to the size of the reduced basis. Hence it already does reflect the 
accuracy of the reduced model. However, the ratio ist not close to one yet, as it would be for a well-conditioned 
surrogate given by the truth-exact evaluation of the residual corresponding to a well-conditioned variational 
formulation. To further improve this ratio by approximating the residual more accurately, we resort to iterative 
tightening as described above. 

The results for a single iteration are recorded in Table [3j It is seen that already after a single run the ratio of 
the surrogate and the true error between the reduced basis approximation and true solution has become much 
closer to one. 

In agreement with the discussion in Section |4.6| the experiments show that a slightly larger number of test 
basis functions than for the convection-diffusion problem is needed here. In particular, unlike in the elliptic 
case non-smooth data (right hand side, boundary shape, and boundary conditions) affect the smoothness of 
the dependence of the solutions on the parameter. In our examples at most a low order polynomial decay of 
the n- widths can be expected. According to Summary |4.11 in Section 4^ since the right hand side is actually 
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dimension 




maximal 


maximal 


error between 


surr / 


trial 


test 


6 


surr 


rb truth 


rb L2 


err 


4 


11 


3.95e-01 


8.44e-03 


2.45e-02 


2.45e-02 


3.45e-01 


6 


17 


4.49e-01 


7.06e-03 


1.40e-02 


1.40e-02 


5.04e-01 


8 


25 


4.87e-01 


4.16e-03 


9.05e-03 


9.05e-03 


4.60e-01 


10 


33 


4.32e-01 


3.37e-03 


5.74e-03 


5.74e-03 


5.87e-01 


12 


40 


4.83e-01 


2.65e-03 


4.65e-03 


4.65e-03 


5.71e-01 


14 


48 


4.23e-01 


1.64e-03 


3.39e-03 


3.39e-03 


4.83e-01 


16 


57 


4.32e-01 


1.50e-03 


2.56e-03 


2.56e-03 


5.84e-01 


18 


65 


4.66e-01 


1.17e-03 


2.33e-03 


2.33e-03 


5.03e-01 


20 


74 


4.16e-01 


1.21e-03 


2.10e-03 


2.10e-03 


5.77e-01 


22 


83 


3.83e-01 


1.02e-03 


1.93e-03 


1.93e-03 


5.29e-01 


24 


91 


4.05e-01 


7.27e-04 


1.58e-03 


1.58e-03 


4.61e-01 



Table 2. Numerical results for the transport problem 6.8 maximal error truth L2 0.000109832 



smooth, we expect that cr„j-, defined in (4.17), that controls the termination of the inner stabilization loop 
drops below the desired 6 < 1 after an acceptable bounded number of steps independent of the dimension of 
the truth space. In fact, one observes that the growth of the test basis stays surprisingly moderate. 





r 



f 



Figure 3. Solution of the transport 



Figure 4. Solution of the transport 



problem (6.8 1, with reduced basis of 
dimension n — 24, m{n) = 91 and an- 
gle n = 0.244579. 



problem (6.9), with reduced basis of 
dimension n = 24, m{n) = 96 and an- 
gle n = 0.256311. 



(6.9) 



Finally, Table |4] and Figure |4] show the results for the problem 
\^sm /iy ^ ^ I X > y 



mn = (0,1) 



P = 



1-y 





X < 0.5 
X > 0.5 



on r_ 



Now the right hand side as well as the boundary conditions exhibit jump discontinuities where the latter is 
transported trough the domain. This causes a further significant reduction of the smoothness of the dependence 



on the solutions on the parameter. Problem (6.9) therefore represents an extreme example involving interacting 



jump discontinuities caused by the right hand side and by the boundary conditions. The small ripples observed 



34 



WOLFGANG DAHMEN, CHRISTIAN PLESKEN, AND GERRIT WELPER 



dimension 
trial 1 test 


S 


maximal 
surr 


maximal 
rb truth 


error between 
rb L2 


surr / 
err 


first reduced basis creation 


20 


81 


3.73e-01 


2.71e-02 


5.46C-02 


5.62e-02 


4.82C-01 


second reduced basis creation 


10 


87 


3.51e-01 


6.45e-02 


7.40e-02 


7.53e-02 


8.57e-01 



Table 3. Numerical results for the transport problem 6.8 after a single cycle of iterative 
tightening. Maximal error truth L2 0.0154814 



dimension 




maximal 


maximal 


error between 


surr / 


trial 


test 


6 


surr 


rb truth 


rb L2 


err 


4 


14 


4.97e-01 


5.91e-02 


1.29e-01 


1.30e-01 


4.54e-01 


6 


23 


4.92e-01 


4.29e-02 


l.OOe-01 


1.02e-01 


4.22e-01 


8 


31 


4.29e-01 


4.34e-02 


7.78e-02 


7.95e-02 


5.46e-01 


10 


40 


4.15e-01 


3.84e-02 


7.78e-02 


7.95e-02 


4.83e-01 


12 


49 


3.71e-01 


3.48e-02 


7.40e-02 


7.53e-02 


4.63e-01 


14 


57 


3.76e-01 


3.12e-02 


6.20e-02 


6.41e-02 


4.87e-01 


16 


64 


3.74e-01 


2.99e-02 


6.20e-02 


6.41e-02 


4.67e-01 


18 


73 


4.63e-01 


2.86e-02 


6.20e-02 


6.41e-02 


4.47e-01 


20 


81 


3.73e-01 


2.71e-02 


5.46e-02 


5.62e-02 


4.82e-01 


22 


87 


4.09e-01 


2.42e-02 


5.46e-02 


5.62e-02 


4.32e-01 


24 


96 


3.91e-01 


2.51e-02 


4.51e-02 


4.79e-02 


5.25e-01 



Table 4. Numerical results for the transport problem 6.9 



maximal error truth L2 0.0154814. 



in the solution plot Figure |4] originate from the superposition of the jumps of the various snapshots involved in 
the solution. Since they do neither grow nor expand one can conclude that the scheme is in fact stable. 

As indicated before, varying the transport direction for such data shows that the dependence of the solution 
on the parameter is even less smooth than in the previous example so that the Kolmogorov widths of the 
solution manifold are expected to decay more slowly. Hence the greedy errors cannot decay too rapidly either. 
Again, by Summary 



4.11 



the quantities tJ„j in (4.17), estimating the number of stabilization steps for X„, 
are expected to decay even more slowly than in the case of zero boundary conditions. Tabel |4] confirms this 
in that slightly more test basis functions are generated than in example |6.8[ Nevertheless, one observes that 
in the initial phase already a few reduced basis functions decrease the error very effectively so that a reduced 
space with trial dimension as low as ten realizes an accuracy that would require a conventional finite element 
space of much larger dimension. Overall, the performance, at least in the given range of truth accuracy, is only 
slightly weaker than for the milder case of zero inflow boundary conditions. The precise implications on the 
approximation of functionals of the solution and possible strategies for alternative ways of enriching the trial 
dictionary will be explored in forthcoming work. 



7. General Saddle Point Problems 

The crucial role of saddle point problems for the generation of well-conditioned variational formulations is 
apparent from the preceding discussion. On the other hand, the concepts developed in this context have an 
immediate bearing on more general saddle point problems of "classical type" . By this we mean (parameter 
dependent) variational formulations e.g. of the Stokes system or those arising in mixed formulations and 
constrained optimization problems. To see this, it is useful to point out the main distinctions between the two 
settings when considering the following general formulation for parameter dependent bilinear forms a^(-,-) : 
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zero boundary 



jump boundary 




- 6-10-2 



4 • 10-2 - 



2 • 10" 



5 10 15 20 

reduced basis trial dimension 




5 10 15 20 25 

reduced basis trial dimension 



Figure 5. Surrogates of the reduced basis approximation for the transport problem (6.8) and (6.9) 



Yf, xYf,^R, bf,{-, ■) : xY^ such that for ^ e P 

-^N afj,{u{iJ,),v) +b{p{^),v) = {f,v), vGY^, 

bi,{q,u{fi)) = {g,q), q€X^. 

For classical problems the following conditions 

(7.2) 

as well as 



(7.3) 
where 



af,{v,w)\<Ca{^i)\\v\\Y^\w\\Y,, \b^.{q,v)\<Cb{^l)\\v\\Y^\q\\x^, v,weY^ q e X^, 
bf,iq,v) 



inf sup 

gex„„gy \\q\\xJv\\Y^ 



> Pifi), a^iv,v) > Cain)\\v\\i, V e V{p), 



Vifi) := {veY^: b^{q, v)^0,yqe M} = kerS* 



are usually met. 



In comparison, the specific structure of the "stabilizing" saddle point problem (3.16) is the following. 



Remark 7.1. is an isomorphism and for a^{v,w) :— {Ry^v,w) one has V{fi) — {0} and condition (7.3) 



holds with Ca(/i) = Ca{^j) — 1 even on Y^. Moreover, when using the norm || • ||^ on X^^, (7.2) holds with 



CbifJ,) — I3{n) =^ 1. However, on the downside, one may encounter failure of either (3.34| or (3.35). 
For the classical problems considered in [131 Ull ESI one can state the following. 



Remark 7.2. The conditions (3.34), (3.35) are both satisfied so that a single reference norm || • ||yxJf can be 



used. Hence, that renormation (3.6) is not necessary for achieving tightness of residual based surrogates which 



now involve both component spaces Y' x X' which can be evaluated by the standard offline- online decomposition, 
see e.g. [Hdg. 



It is well known (see e.g. |2) that, given (7.2 1 and (7.3), the validity of the mapping property MP and the 
best approximation property BAP hinges again on the inf-sup condition 



(7.4) 



^Ht '^^P Til — iTii — - py-.wW- 



36 



WOLFGANG DAHMEN, CHRISTIAN PLESKEN, AND GERRIT WELPER 



For classical problems it is known for V = Y,W — X, has to be ensured for the truth spaces V — Yj^, W — Xj^ 
through suitably chosen finite element spaces, say, and again need to be ensured by stabilizing strategies for the 
reduced spaces V = Yn, W — Xn- 

In view of Remark [7.2[ both schemes Update-(5 and Update-inf-SUP can be applied. Since the spaces F„ no 
longer just serve as stabilizers but need to contribute to the target approximation accuracy of the full solution 
manifold 



(7.5) M := {[u{n),p{n)] : solves JtI]), fi e P} =: My x Mx, 

the only changes that need to be incorporated in a slightly modified version DG-2 of DG-1 are: 

• In Algorithm Q replace step 5 by: 
Set 

span {Xn ,p} ^ Xn, span {y„ , u} F„ , 

i.e., both component spaces are updated in the outer greedy step. 

• In Algorithm ([5| step 5 is replaced by 
Yn.Xn ^ Update- AppROxiMATiON(y„, X„) 

• Replace the surrogate by 

R*{^i,V->< W) := \\f ~ A^uw,v{^l) - B^,pwy{tJ^)\\Y;^ + h - B*uw,v{p)\\x'^, 
see [HITl]. 



Clearly, under the given assumptions (3.34), (3.351, is compact. Denoting again by pn{p-) , Un{p-) the 



solution components produced by the scheme DG-2 and comparing the greedy errors (5.171 



<7n{M) := sup {IIp(m) ~'Pni^J■)\\x + - Un{p)\\Y}, 

i-iev 



with the n-widths dn{M)xy.Y and keeping Proposition 4.10 in mind, we extend the results in [131 [111 [26] as 
follows. 



Corollary 7.3. The scheme DG-2 applied to (7.1) is under the above assumptions rate-optimal. 



8. Concluding Remarks 

The generation of well-conditioned variational formulations for non-coercive or indefinite problems has been 
proposed as the central ingredient of a general strategy for constructing tight surrogates for RBMs also for 
such problem classes. In contrast to previous work, well-conditioned tight surrogates are obtained in a feasible 
way in all settings warranting a near-optimal performance of the corresponding REM, which does not seem to 
be achievable with the aid of previously known concepts. We emphasize that these concepts apply as well to 
space-time discretizations of unsteady problems (see [7] ) offering interesting perspectives with regard to robustly 
capturing long-term dynamics. The presented application to two simple model problems is to be viewed as a 
first proof of concept. The two examples are to bring out some essential obstructions and raise issues that have 
so far not been addressed in this context. In particular, they hint at the principal limitations of RBMs in their 
standard formulations, especially regarding the smoothness of the parameter dependence. 
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